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 #include "G4HepRepFileSceneHandler.hh"
00034 #include "G4HepRepFile.hh"
00035 #include "G4HepRepMessenger.hh"
00036 #include "G4UIcommand.hh"
00037
00038 #include "G4PhysicalConstants.hh"
00039 #include "G4SystemOfUnits.hh"
00040 #include "G4Version.hh"
00041 #include "G4VSolid.hh"
00042 #include "G4PhysicalVolumeModel.hh"
00043 #include "G4VPhysicalVolume.hh"
00044 #include "G4LogicalVolume.hh"
00045 #include "G4RotationMatrix.hh"
00046 #include "G4Material.hh"
00047 #include "G4Polymarker.hh"
00048 #include "G4Polyline.hh"
00049 #include "G4Text.hh"
00050 #include "G4Circle.hh"
00051 #include "G4Square.hh"
00052 #include "G4Polyhedron.hh"
00053 #include "G4NURBS.hh"
00054 #include "G4VTrajectory.hh"
00055 #include "G4VTrajectoryPoint.hh"
00056 #include "G4TrajectoriesModel.hh"
00057 #include "G4VHit.hh"
00058 #include "G4AttDef.hh"
00059 #include "G4AttValue.hh"
00060 #include "G4AttCheck.hh"
00061 #include "G4VisManager.hh"
00062 #include "G4VisTrajContext.hh"
00063 #include "G4VTrajectoryModel.hh"
00064
00065
00066 #include "G4HepRepFileXMLWriter.hh"
00067
00068 G4int G4HepRepFileSceneHandler::fSceneIdCount = 0;
00069
00070
00071 G4HepRepFileSceneHandler::G4HepRepFileSceneHandler(G4VGraphicsSystem& system,
00072 const G4String& name):
00073 G4VSceneHandler(system, fSceneIdCount++, name)
00074 {
00075 hepRepXMLWriter = ((G4HepRepFile*)(&system))->GetHepRepXMLWriter();
00076 fileCounter = 0;
00077
00078 inPrimitives2D = false;
00079 warnedAbout3DText = false;
00080 warnedAbout2DMarkers = false;
00081 haveVisible = false;
00082 drawingTraj = false;
00083 doneInitTraj = false;
00084 drawingHit = false;
00085 doneInitHit = false;
00086 }
00087
00088
00089 G4HepRepFileSceneHandler::~G4HepRepFileSceneHandler() {}
00090
00091
00092 void G4HepRepFileSceneHandler::BeginModeling() {
00093 G4VisManager* visManager = G4VisManager::GetInstance();
00094 const G4VTrajectoryModel* model = visManager->CurrentTrajDrawModel();
00095 trajContext = & model->GetContext();
00096
00097 G4VSceneHandler::BeginModeling();
00098 }
00099
00100
00101 void G4HepRepFileSceneHandler::EndModeling() {
00102 G4VSceneHandler::EndModeling();
00103 }
00104
00105 void G4HepRepFileSceneHandler::BeginPrimitives2D
00106 (const G4Transform3D& objectTransformation) {
00107 #ifdef G4HEPREPFILEDEBUG
00108 G4cout << "G4HepRepFileSceneHandler::BeginPrimitives2D() " << G4endl;
00109 #endif
00110 inPrimitives2D = true;
00111 G4VSceneHandler::BeginPrimitives2D(objectTransformation);
00112 }
00113
00114 void G4HepRepFileSceneHandler::EndPrimitives2D () {
00115 #ifdef G4HEPREPFILEDEBUG
00116 G4cout << "G4HepRepFileSceneHandler::EndPrimitives2D() " << G4endl;
00117 #endif
00118 G4VSceneHandler::EndPrimitives2D();
00119 inPrimitives2D = false;
00120 }
00121
00122
00123 #ifdef G4HEPREPFILEDEBUG
00124 void G4HepRepFileSceneHandler::PrintThings() {
00125 G4cout <<
00126 " with transformation "
00127 << fObjectTransformation;
00128 G4PhysicalVolumeModel* pPVModel =
00129 dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
00130 if (pPVModel) {
00131 G4VPhysicalVolume* pCurrentPV = pPVModel->GetCurrentPV();
00132 G4LogicalVolume* pCurrentLV = pPVModel->GetCurrentLV();
00133 G4int currentDepth = pPVModel->GetCurrentDepth();
00134 G4cout <<
00135 "\n current physical volume: "
00136 << pCurrentPV->GetName() <<
00137 "\n current logical volume: "
00138 << pCurrentLV->GetName() <<
00139 "\n current depth of geometry tree: "
00140 << currentDepth;
00141 }
00142 G4cout << G4endl;
00143 }
00144 #endif
00145
00146
00147 void G4HepRepFileSceneHandler::AddSolid(const G4Box& box) {
00148 #ifdef G4HEPREPFILEDEBUG
00149 G4cout <<
00150 "G4HepRepFileSceneHandler::AddSolid(const G4Box& box) called for "
00151 << box.GetName()
00152 << G4endl;
00153 PrintThings();
00154 #endif
00155
00156 if (drawingTraj)
00157 return;
00158
00159 if (drawingHit)
00160 InitHit();
00161
00162 haveVisible = false;
00163 AddHepRepInstance("Prism", NULL);
00164
00165 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00166
00167 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
00168 return;
00169
00170 hepRepXMLWriter->addPrimitive();
00171
00172 G4double dx = box.GetXHalfLength();
00173 G4double dy = box.GetYHalfLength();
00174 G4double dz = box.GetZHalfLength();
00175
00176 G4Point3D vertex1(G4Point3D( dx, dy,-dz));
00177 G4Point3D vertex2(G4Point3D( dx,-dy,-dz));
00178 G4Point3D vertex3(G4Point3D(-dx,-dy,-dz));
00179 G4Point3D vertex4(G4Point3D(-dx, dy,-dz));
00180 G4Point3D vertex5(G4Point3D( dx, dy, dz));
00181 G4Point3D vertex6(G4Point3D( dx,-dy, dz));
00182 G4Point3D vertex7(G4Point3D(-dx,-dy, dz));
00183 G4Point3D vertex8(G4Point3D(-dx, dy, dz));
00184
00185 vertex1 = (fObjectTransformation) * vertex1;
00186 vertex2 = (fObjectTransformation) * vertex2;
00187 vertex3 = (fObjectTransformation) * vertex3;
00188 vertex4 = (fObjectTransformation) * vertex4;
00189 vertex5 = (fObjectTransformation) * vertex5;
00190 vertex6 = (fObjectTransformation) * vertex6;
00191 vertex7 = (fObjectTransformation) * vertex7;
00192 vertex8 = (fObjectTransformation) * vertex8;
00193
00194 hepRepXMLWriter->addPoint(vertex1.x(), vertex1.y(), vertex1.z());
00195 hepRepXMLWriter->addPoint(vertex2.x(), vertex2.y(), vertex2.z());
00196 hepRepXMLWriter->addPoint(vertex3.x(), vertex3.y(), vertex3.z());
00197 hepRepXMLWriter->addPoint(vertex4.x(), vertex4.y(), vertex4.z());
00198 hepRepXMLWriter->addPoint(vertex5.x(), vertex5.y(), vertex5.z());
00199 hepRepXMLWriter->addPoint(vertex6.x(), vertex6.y(), vertex6.z());
00200 hepRepXMLWriter->addPoint(vertex7.x(), vertex7.y(), vertex7.z());
00201 hepRepXMLWriter->addPoint(vertex8.x(), vertex8.y(), vertex8.z());
00202 }
00203
00204
00205 void G4HepRepFileSceneHandler::AddSolid(const G4Cons& cons) {
00206 #ifdef G4HEPREPFILEDEBUG
00207 G4cout <<
00208 "G4HepRepFileSceneHandler::AddSolid(const G4Cons& cons) called for "
00209 << cons.GetName()
00210 << G4endl;
00211 PrintThings();
00212 #endif
00213
00214
00215
00216 G4RotationMatrix r = fObjectTransformation.getRotation();
00217 G4bool linedUpWithAnAxis = (std::fabs(r.phiX())<=.001 ||
00218 std::fabs(r.phiY())<=.001 ||
00219 std::fabs(r.phiZ())<=.001 ||
00220 std::fabs(r.phiX()-pi)<=.001 ||
00221 std::fabs(r.phiY()-pi)<=.001 ||
00222 std::fabs(r.phiZ()-pi)<=.001);
00223
00224
00225
00226
00227
00228
00229 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00230 if (cons.GetDeltaPhiAngle() < twopi || !linedUpWithAnAxis || messenger->renderCylAsPolygons())
00231 {
00232 G4VSceneHandler::AddSolid(cons);
00233 } else {
00234
00235 if (drawingTraj)
00236 return;
00237
00238 if (drawingHit)
00239 InitHit();
00240
00241 haveVisible = false;
00242 AddHepRepInstance("Cylinder", NULL);
00243
00244 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
00245 return;
00246
00247 G4Point3D vertex1(G4Point3D( 0., 0., -cons.GetZHalfLength()));
00248 G4Point3D vertex2(G4Point3D( 0., 0., cons.GetZHalfLength()));
00249
00250 vertex1 = (fObjectTransformation) * vertex1;
00251 vertex2 = (fObjectTransformation) * vertex2;
00252
00253
00254 hepRepXMLWriter->addPrimitive();
00255 hepRepXMLWriter->addAttValue("Radius1",messenger->getScale() * cons.GetOuterRadiusMinusZ());
00256 hepRepXMLWriter->addAttValue("Radius2",messenger->getScale() * cons.GetOuterRadiusPlusZ());
00257 hepRepXMLWriter->addPoint(vertex1.x(), vertex1.y(), vertex1.z());
00258 hepRepXMLWriter->addPoint(vertex2.x(), vertex2.y(), vertex2.z());
00259
00260
00261 hepRepXMLWriter->addPrimitive();
00262 hepRepXMLWriter->addAttValue("Radius1",messenger->getScale() * cons.GetInnerRadiusMinusZ());
00263 hepRepXMLWriter->addAttValue("Radius2",messenger->getScale() * cons.GetInnerRadiusPlusZ());
00264 hepRepXMLWriter->addPoint(vertex1.x(), vertex1.y(), vertex1.z());
00265 hepRepXMLWriter->addPoint(vertex2.x(), vertex2.y(), vertex2.z());
00266 }
00267 }
00268
00269
00270 void G4HepRepFileSceneHandler::AddSolid(const G4Tubs& tubs) {
00271 #ifdef G4HEPREPFILEDEBUG
00272 G4cout <<
00273 "G4HepRepFileSceneHandler::AddSolid(const G4Tubs& tubs) called for "
00274 << tubs.GetName()
00275 << G4endl;
00276 PrintThings();
00277 #endif
00278
00279
00280
00281 G4RotationMatrix r = fObjectTransformation.getRotation();
00282 G4bool linedUpWithAnAxis = (std::fabs(r.phiX())<=.001 ||
00283 std::fabs(r.phiY())<=.001 ||
00284 std::fabs(r.phiZ())<=.001 ||
00285 std::fabs(r.phiX()-pi)<=.001 ||
00286 std::fabs(r.phiY()-pi)<=.001 ||
00287 std::fabs(r.phiZ()-pi)<=.001);
00288
00289
00290
00291
00292
00293
00294 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00295 if (tubs.GetDeltaPhiAngle() < twopi || !linedUpWithAnAxis || messenger->renderCylAsPolygons())
00296 {
00297 G4VSceneHandler::AddSolid(tubs);
00298 } else {
00299
00300 if (drawingTraj)
00301 return;
00302
00303 if (drawingHit)
00304 InitHit();
00305
00306 haveVisible = false;
00307 AddHepRepInstance("Cylinder", NULL);
00308
00309 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
00310 return;
00311
00312 G4Point3D vertex1(G4Point3D( 0., 0., -tubs.GetZHalfLength()));
00313 G4Point3D vertex2(G4Point3D( 0., 0., tubs.GetZHalfLength()));
00314
00315 vertex1 = (fObjectTransformation) * vertex1;
00316 vertex2 = (fObjectTransformation) * vertex2;
00317
00318
00319 hepRepXMLWriter->addPrimitive();
00320 hepRepXMLWriter->addAttValue("Radius1", messenger->getScale() * tubs.GetOuterRadius());
00321 hepRepXMLWriter->addAttValue("Radius2", messenger->getScale() * tubs.GetOuterRadius());
00322 hepRepXMLWriter->addPoint(vertex1.x(), vertex1.y(), vertex1.z());
00323 hepRepXMLWriter->addPoint(vertex2.x(), vertex2.y(), vertex2.z());
00324
00325
00326 if (tubs.GetInnerRadius() != 0.) {
00327 hepRepXMLWriter->addPrimitive();
00328 hepRepXMLWriter->addAttValue("Radius1", messenger->getScale() * tubs.GetInnerRadius());
00329 hepRepXMLWriter->addAttValue("Radius2", messenger->getScale() * tubs.GetInnerRadius());
00330 hepRepXMLWriter->addPoint(vertex1.x(), vertex1.y(), vertex1.z());
00331 hepRepXMLWriter->addPoint(vertex2.x(), vertex2.y(), vertex2.z());
00332 }
00333 }
00334 }
00335
00336
00337 void G4HepRepFileSceneHandler::AddSolid(const G4Trd& trd) {
00338 #ifdef G4HEPREPFILEDEBUG
00339 G4cout <<
00340 "G4HepRepFileSceneHandler::AddSolid(const G4Trd& trd) called for "
00341 << trd.GetName()
00342 << G4endl;
00343 PrintThings();
00344 #endif
00345
00346 if (drawingTraj)
00347 return;
00348
00349 if (drawingHit)
00350 InitHit();
00351
00352 haveVisible = false;
00353 AddHepRepInstance("Prism", NULL);
00354
00355 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00356
00357 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
00358 return;
00359
00360 hepRepXMLWriter->addPrimitive();
00361
00362 G4double dx1 = trd.GetXHalfLength1();
00363 G4double dy1 = trd.GetYHalfLength1();
00364 G4double dx2 = trd.GetXHalfLength2();
00365 G4double dy2 = trd.GetYHalfLength2();
00366 G4double dz = trd.GetZHalfLength();
00367
00368 G4Point3D vertex1(G4Point3D( dx1, dy1,-dz));
00369 G4Point3D vertex2(G4Point3D( dx1,-dy1,-dz));
00370 G4Point3D vertex3(G4Point3D(-dx1,-dy1,-dz));
00371 G4Point3D vertex4(G4Point3D(-dx1, dy1,-dz));
00372 G4Point3D vertex5(G4Point3D( dx2, dy2, dz));
00373 G4Point3D vertex6(G4Point3D( dx2,-dy2, dz));
00374 G4Point3D vertex7(G4Point3D(-dx2,-dy2, dz));
00375 G4Point3D vertex8(G4Point3D(-dx2, dy2, dz));
00376
00377 vertex1 = (fObjectTransformation) * vertex1;
00378 vertex2 = (fObjectTransformation) * vertex2;
00379 vertex3 = (fObjectTransformation) * vertex3;
00380 vertex4 = (fObjectTransformation) * vertex4;
00381 vertex5 = (fObjectTransformation) * vertex5;
00382 vertex6 = (fObjectTransformation) * vertex6;
00383 vertex7 = (fObjectTransformation) * vertex7;
00384 vertex8 = (fObjectTransformation) * vertex8;
00385
00386 hepRepXMLWriter->addPoint(vertex1.x(), vertex1.y(), vertex1.z());
00387 hepRepXMLWriter->addPoint(vertex2.x(), vertex2.y(), vertex2.z());
00388 hepRepXMLWriter->addPoint(vertex3.x(), vertex3.y(), vertex3.z());
00389 hepRepXMLWriter->addPoint(vertex4.x(), vertex4.y(), vertex4.z());
00390 hepRepXMLWriter->addPoint(vertex5.x(), vertex5.y(), vertex5.z());
00391 hepRepXMLWriter->addPoint(vertex6.x(), vertex6.y(), vertex6.z());
00392 hepRepXMLWriter->addPoint(vertex7.x(), vertex7.y(), vertex7.z());
00393 hepRepXMLWriter->addPoint(vertex8.x(), vertex8.y(), vertex8.z());
00394 }
00395
00396
00397 void G4HepRepFileSceneHandler::AddSolid(const G4Trap& trap) {
00398 #ifdef G4HEPREPFILEDEBUG
00399 G4cout <<
00400 "G4HepRepFileSceneHandler::AddSolid(const G4Trap& trap) called for "
00401 << trap.GetName()
00402 << G4endl;
00403 PrintThings();
00404 #endif
00405 G4VSceneHandler::AddSolid(trap);
00406 }
00407
00408
00409 void G4HepRepFileSceneHandler::AddSolid(const G4Sphere& sphere) {
00410 #ifdef G4HEPREPFILEDEBUG
00411 G4cout <<
00412 "G4HepRepFileSceneHandler::AddSolid(const G4Sphere& sphere) called for "
00413 << sphere.GetName()
00414 << G4endl;
00415 PrintThings();
00416 #endif
00417 G4VSceneHandler::AddSolid(sphere);
00418 }
00419
00420
00421 void G4HepRepFileSceneHandler::AddSolid(const G4Para& para) {
00422 #ifdef G4HEPREPFILEDEBUG
00423 G4cout <<
00424 "G4HepRepFileSceneHandler::AddSolid(const G4Para& para) called for "
00425 << para.GetName()
00426 << G4endl;
00427 PrintThings();
00428 #endif
00429 G4VSceneHandler::AddSolid(para);
00430 }
00431
00432
00433 void G4HepRepFileSceneHandler::AddSolid(const G4Torus& torus) {
00434 #ifdef G4HEPREPFILEDEBUG
00435 G4cout <<
00436 "G4HepRepFileSceneHandler::AddSolid(const G4Torus& torus) called for "
00437 << torus.GetName()
00438 << G4endl;
00439 PrintThings();
00440 #endif
00441 G4VSceneHandler::AddSolid(torus);
00442 }
00443
00444
00445 void G4HepRepFileSceneHandler::AddSolid(const G4Polycone& polycone) {
00446 #ifdef G4HEPREPFILEDEBUG
00447 G4cout <<
00448 "G4HepRepFileSceneHandler::AddSolid(const G4Polycone& polycone) called for "
00449 << polycone.GetName()
00450 << G4endl;
00451 PrintThings();
00452 #endif
00453 G4VSceneHandler::AddSolid(polycone);
00454 }
00455
00456
00457 void G4HepRepFileSceneHandler::AddSolid(const G4Polyhedra& polyhedra) {
00458 #ifdef G4HEPREPFILEDEBUG
00459 G4cout <<
00460 "G4HepRepFileSceneHandler::AddSolid(const G4Polyhedra& polyhedra) called for "
00461 << polyhedra.GetName()
00462 << G4endl;
00463 PrintThings();
00464 #endif
00465 G4VSceneHandler::AddSolid(polyhedra);
00466 }
00467
00468
00469 void G4HepRepFileSceneHandler::AddSolid(const G4VSolid& solid) {
00470 #ifdef G4HEPREPFILEDEBUG
00471 G4cout <<
00472 "G4HepRepFileSceneHandler::AddSolid(const G4Solid& solid) called for "
00473 << solid.GetName()
00474 << G4endl;
00475 PrintThings();
00476 #endif
00477 G4VSceneHandler::AddSolid(solid);
00478 }
00479
00480
00481 void G4HepRepFileSceneHandler::AddCompound (const G4VTrajectory& traj) {
00482 #ifdef G4HEPREPFILEDEBUG
00483 G4cout << "G4HepRepFileSceneHandler::AddCompound(const G4VTrajectory&) " << G4endl;
00484 #endif
00485
00486 G4TrajectoriesModel* pTrModel =
00487 dynamic_cast<G4TrajectoriesModel*>(fpModel);
00488 if (!pTrModel) G4Exception
00489 ("G4HepRepFileSceneHandler::AddCompound(const G4VTrajectory&)",
00490 "vis-HepRep0001", FatalException, "Not a G4TrajectoriesModel.");
00491
00492
00493 std::vector<G4AttValue>* rawTrajAttValues = traj.CreateAttValues();
00494 trajAttValues =
00495 new std::vector<G4AttValue>;
00496 trajAttDefs =
00497 new std::map<G4String,G4AttDef>;
00498
00499
00500 std::vector<G4AttValue>::iterator iAttVal;
00501 std::map<G4String,G4AttDef>::const_iterator iAttDef;
00502 G4int i;
00503
00504
00505
00506 if (rawTrajAttValues) {
00507 G4bool error = G4AttCheck(rawTrajAttValues,
00508 traj.GetAttDefs()).Standard(trajAttValues,trajAttDefs);
00509 if (error) {
00510 G4cout << "G4HepRepFileSceneHandler::AddCompound(traj):"
00511 "\nERROR found during conversion to standard trajectory attributes."
00512 << G4endl;
00513 }
00514 #ifdef G4HEPREPFILEDEBUG
00515 G4cout <<
00516 "G4HepRepFileSceneHandler::AddCompound(traj): standardised attributes:\n"
00517 << G4AttCheck(trajAttValues,trajAttDefs) << G4endl;
00518 #endif
00519 delete rawTrajAttValues;
00520 }
00521
00522
00523 CheckFileOpen();
00524
00525
00526 if (strcmp("Event Data",hepRepXMLWriter->prevTypeName[0])!=0) {
00527 hepRepXMLWriter->addType("Event Data",0);
00528 hepRepXMLWriter->addInstance();
00529 }
00530
00531
00532 G4String previousName = hepRepXMLWriter->prevTypeName[1];
00533 hepRepXMLWriter->addType("Trajectories",1);
00534
00535
00536
00537 if (strcmp("Trajectories",previousName)!=0) {
00538 hepRepXMLWriter->addAttValue("Layer",100);
00539
00540
00541
00542
00543
00544 if (trajAttValues && trajAttDefs) {
00545 for (iAttVal = trajAttValues->begin();
00546 iAttVal != trajAttValues->end(); ++iAttVal) {
00547 iAttDef = trajAttDefs->find(iAttVal->GetName());
00548 if (iAttDef != trajAttDefs->end()) {
00549
00550
00551 G4String category = iAttDef->second.GetCategory();
00552 if (strcmp(category,"Draw")!=0 &&
00553 strcmp(category,"Physics")!=0 &&
00554 strcmp(category,"Association")!=0 &&
00555 strcmp(category,"PickAction")!=0)
00556 category = "Physics";
00557 hepRepXMLWriter->addAttDef(iAttVal->GetName(), iAttDef->second.GetDesc(),
00558 category, iAttDef->second.GetExtra());
00559 }
00560 }
00561 }
00562
00563
00564
00565
00566 if ((trajContext->GetDrawStepPts() || trajContext->GetDrawAuxPts())
00567 && traj.GetPointEntries()>0) {
00568 G4VTrajectoryPoint* aTrajectoryPoint = traj.GetPoint(0);
00569
00570
00571 std::vector<G4AttValue>* rawPointAttValues = aTrajectoryPoint->CreateAttValues();
00572 std::vector<G4AttValue>* pointAttValues =
00573 new std::vector<G4AttValue>;
00574 std::map<G4String,G4AttDef>* pointAttDefs =
00575 new std::map<G4String,G4AttDef>;
00576
00577
00578
00579 if (rawPointAttValues) {
00580 G4bool error = G4AttCheck(rawPointAttValues,
00581 aTrajectoryPoint->GetAttDefs()).Standard(pointAttValues,pointAttDefs);
00582 if (error) {
00583 G4cout << "G4HepRepFileSceneHandler::AddCompound(traj):"
00584 "\nERROR found during conversion to standard first point attributes." << G4endl;
00585 }
00586
00587
00588 if (pointAttValues && pointAttDefs) {
00589 for (iAttVal = pointAttValues->begin();
00590 iAttVal != pointAttValues->end(); ++iAttVal) {
00591 iAttDef =
00592 pointAttDefs->find(iAttVal->GetName());
00593 if (iAttDef != pointAttDefs->end()) {
00594
00595
00596 G4String category = iAttDef->second.GetCategory();
00597 if (strcmp(category,"Draw")!=0 &&
00598 strcmp(category,"Physics")!=0 &&
00599 strcmp(category,"Association")!=0 &&
00600 strcmp(category,"PickAction")!=0)
00601 category = "Physics";
00602
00603
00604
00605 if (strcmp(iAttVal->GetName(),"Aux-X")!=0 &&
00606 strcmp(iAttVal->GetName(),"Aux-Y")!=0 &&
00607 strcmp(iAttVal->GetName(),"Aux-Z")!=0 &&
00608 strcmp(iAttVal->GetName(),"Pos-X")!=0 &&
00609 strcmp(iAttVal->GetName(),"Pos-Y")!=0 &&
00610 strcmp(iAttVal->GetName(),"Pos-Z")!=0)
00611 hepRepXMLWriter->addAttDef(iAttVal->GetName(), iAttDef->second.GetDesc(),
00612 category, iAttDef->second.GetExtra());
00613 }
00614 }
00615 }
00616 delete rawPointAttValues;
00617 }
00618
00619
00620 if (pointAttValues)
00621 delete pointAttValues;
00622 if (pointAttDefs)
00623 delete pointAttDefs;
00624 }
00625 }
00626
00627
00628
00629
00630
00631
00632
00633
00634
00635 drawingTraj = true;
00636 doneInitTraj = false;
00637 G4VSceneHandler::AddCompound(traj);
00638 drawingTraj = false;
00639
00640
00641 if (trajContext->GetDrawStepPts()) {
00642 if (!doneInitTraj)
00643 InitTrajectory();
00644
00645
00646
00647
00648 previousName = hepRepXMLWriter->prevTypeName[2];
00649 hepRepXMLWriter->addType("Trajectory Step Points",2);
00650
00651 float redness;
00652 float greenness;
00653 float blueness;
00654 G4int markSize;
00655 G4bool visible;
00656 G4bool square;
00657 G4Colour colour = trajContext->GetStepPtsColour();
00658 redness = colour.GetRed();
00659 greenness = colour.GetGreen();
00660 blueness = colour.GetBlue();
00661 markSize = (G4int) trajContext->GetStepPtsSize();
00662 visible = (G4int) trajContext->GetStepPtsVisible();
00663 square = (trajContext->GetStepPtsType()==G4Polymarker::squares);
00664
00665
00666 if (redness==0. && greenness==0. && blueness==0.) {
00667 redness = 1.;
00668 greenness = 1.;
00669 blueness = 1.;
00670 }
00671
00672
00673 if (strcmp("Trajectory Step Points",previousName)!=0) {
00674 hepRepXMLWriter->addAttValue("DrawAs","Point");
00675 hepRepXMLWriter->addAttValue("MarkColor", redness, greenness, blueness);
00676 hepRepXMLWriter->addAttValue("MarkSize",markSize);
00677 hepRepXMLWriter->addAttValue("Layer",110);
00678 hepRepXMLWriter->addAttValue("Visibility",visible);
00679 if (square)
00680 hepRepXMLWriter->addAttValue("MarkName","square");
00681 else
00682 hepRepXMLWriter->addAttValue("MarkName","dot");
00683 }
00684
00685
00686 for (i = 0; i < traj.GetPointEntries(); i++) {
00687 G4VTrajectoryPoint* aTrajectoryPoint = traj.GetPoint(i);
00688
00689
00690 hepRepXMLWriter->addInstance();
00691
00692
00693 std::vector<G4AttValue>* rawPointAttValues = aTrajectoryPoint->CreateAttValues();
00694 std::vector<G4AttValue>* pointAttValues =
00695 new std::vector<G4AttValue>;
00696 std::map<G4String,G4AttDef>* pointAttDefs =
00697 new std::map<G4String,G4AttDef>;
00698
00699
00700
00701 if (rawPointAttValues) {
00702 G4bool error = G4AttCheck(rawPointAttValues,
00703 aTrajectoryPoint->GetAttDefs()).Standard(pointAttValues,pointAttDefs);
00704 if (error) {
00705 G4cout << "G4HepRepFileSceneHandler::AddCompound(traj):"
00706 "\nERROR found during conversion to standard point attributes." << G4endl;
00707 }
00708
00709
00710 if (pointAttValues) {
00711 for (iAttVal = pointAttValues->begin();
00712 iAttVal != pointAttValues->end(); ++iAttVal)
00713
00714
00715
00716 if (strcmp(iAttVal->GetName(),"Aux-X")!=0 &&
00717 strcmp(iAttVal->GetName(),"Aux-Y")!=0 &&
00718 strcmp(iAttVal->GetName(),"Aux-Z")!=0 &&
00719 strcmp(iAttVal->GetName(),"Pos-X")!=0 &&
00720 strcmp(iAttVal->GetName(),"Pos-Y")!=0 &&
00721 strcmp(iAttVal->GetName(),"Pos-Z")!=0)
00722 hepRepXMLWriter->addAttValue(iAttVal->GetName(), iAttVal->GetValue());
00723 delete pointAttValues;
00724 }
00725 delete rawPointAttValues;
00726 }
00727
00728
00729 if (pointAttDefs)
00730 delete pointAttDefs;
00731
00732
00733 hepRepXMLWriter->addPrimitive();
00734 G4Point3D vertex = aTrajectoryPoint->GetPosition();
00735 hepRepXMLWriter->addPoint(vertex.x(), vertex.y(), vertex.z());
00736 }
00737 }
00738
00739
00740 if (trajContext->GetDrawAuxPts()) {
00741 if (!doneInitTraj)
00742 InitTrajectory();
00743
00744
00745
00746
00747 previousName = hepRepXMLWriter->prevTypeName[2];
00748 hepRepXMLWriter->addType("Trajectory Auxiliary Points",2);
00749
00750 float redness;
00751 float greenness;
00752 float blueness;
00753 G4int markSize;
00754 G4bool visible;
00755 G4bool square;
00756 G4Colour colour = trajContext->GetAuxPtsColour();
00757 redness = colour.GetRed();
00758 greenness = colour.GetGreen();
00759 blueness = colour.GetBlue();
00760 markSize = (G4int) trajContext->GetAuxPtsSize();
00761 visible = (G4int) trajContext->GetAuxPtsVisible();
00762 square = (trajContext->GetAuxPtsType()==G4Polymarker::squares);
00763
00764
00765 if (redness==0. && greenness==0. && blueness==0.) {
00766 redness = 1.;
00767 greenness = 1.;
00768 blueness = 1.;
00769 }
00770
00771
00772 if (strcmp("Trajectory Auxiliary Points",previousName)!=0) {
00773 hepRepXMLWriter->addAttValue("DrawAs","Point");
00774 hepRepXMLWriter->addAttValue("MarkColor", redness, greenness, blueness);
00775 hepRepXMLWriter->addAttValue("MarkSize",markSize);
00776 hepRepXMLWriter->addAttValue("Layer",110);
00777 hepRepXMLWriter->addAttValue("Visibility",visible);
00778 if (square)
00779 hepRepXMLWriter->addAttValue("MarkName","Square");
00780 else
00781 hepRepXMLWriter->addAttValue("MarkName","Dot");
00782 }
00783
00784
00785 for (i = 0; i < traj.GetPointEntries(); i++) {
00786 G4VTrajectoryPoint* aTrajectoryPoint = traj.GetPoint(i);
00787
00788
00789 hepRepXMLWriter->addInstance();
00790
00791
00792 std::vector<G4AttValue>* rawPointAttValues = aTrajectoryPoint->CreateAttValues();
00793 std::vector<G4AttValue>* pointAttValues =
00794 new std::vector<G4AttValue>;
00795 std::map<G4String,G4AttDef>* pointAttDefs =
00796 new std::map<G4String,G4AttDef>;
00797
00798
00799
00800 if (rawPointAttValues) {
00801 G4bool error = G4AttCheck(rawPointAttValues,
00802 aTrajectoryPoint->GetAttDefs()).Standard(pointAttValues,pointAttDefs);
00803 if (error) {
00804 G4cout << "G4HepRepFileSceneHandler::AddCompound(traj):"
00805 "\nERROR found during conversion to standard point attributes." << G4endl;
00806 }
00807
00808
00809 if (pointAttValues) {
00810 for (iAttVal = pointAttValues->begin();
00811 iAttVal != pointAttValues->end(); ++iAttVal)
00812
00813
00814
00815 if (strcmp(iAttVal->GetName(),"Aux-X")!=0 &&
00816 strcmp(iAttVal->GetName(),"Aux-Y")!=0 &&
00817 strcmp(iAttVal->GetName(),"Aux-Z")!=0 &&
00818 strcmp(iAttVal->GetName(),"Pos-X")!=0 &&
00819 strcmp(iAttVal->GetName(),"Pos-Y")!=0 &&
00820 strcmp(iAttVal->GetName(),"Pos-Z")!=0)
00821 hepRepXMLWriter->addAttValue(iAttVal->GetName(), iAttVal->GetValue());
00822 delete pointAttValues;
00823 }
00824 delete rawPointAttValues;
00825 }
00826
00827
00828 if (pointAttDefs)
00829 delete pointAttDefs;
00830
00831
00832 G4Point3D vertex = aTrajectoryPoint->GetPosition();
00833
00834
00835 const std::vector<G4ThreeVector>* auxiliaries = aTrajectoryPoint->GetAuxiliaryPoints();
00836 if (0 != auxiliaries) {
00837 for (size_t iAux=0; iAux<auxiliaries->size(); ++iAux) {
00838 const G4ThreeVector auxPos((*auxiliaries)[iAux]);
00839 hepRepXMLWriter->addPrimitive();
00840 hepRepXMLWriter->addPoint(auxPos.x(), auxPos.y(), auxPos.z());
00841 }
00842 }
00843 }
00844 }
00845 }
00846
00847
00848 void G4HepRepFileSceneHandler::AddCompound (const G4VHit& hit) {
00849 #ifdef G4HEPREPFILEDEBUG
00850 G4cout << "G4HepRepFileSceneHandler::AddCompound(G4VHit&) " << G4endl;
00851 #endif
00852
00853
00854 std::vector<G4AttValue>* rawHitAttValues = hit.CreateAttValues();
00855 hitAttValues =
00856 new std::vector<G4AttValue>;
00857 hitAttDefs =
00858 new std::map<G4String,G4AttDef>;
00859
00860
00861 std::vector<G4AttValue>::iterator iAttVal;
00862 std::map<G4String,G4AttDef>::const_iterator iAttDef;
00863
00864
00865
00866 if (rawHitAttValues) {
00867 G4bool error = G4AttCheck(rawHitAttValues,
00868 hit.GetAttDefs()).Standard(hitAttValues,hitAttDefs);
00869 if (error) {
00870 G4cout << "G4HepRepFileSceneHandler::AddCompound(hit):"
00871 "\nERROR found during conversion to standard hit attributes."
00872 << G4endl;
00873 }
00874 #ifdef G4HEPREPFILEDEBUG
00875 G4cout <<
00876 "G4HepRepFileSceneHandler::AddCompound(hit): standardised attributes:\n"
00877 << G4AttCheck(hitAttValues,hitAttDefs) << G4endl;
00878 #endif
00879 delete rawHitAttValues;
00880 }
00881
00882
00883 CheckFileOpen();
00884
00885
00886 if (strcmp("Event Data",hepRepXMLWriter->prevTypeName[0])!=0) {
00887 hepRepXMLWriter->addType("Event Data",0);
00888 hepRepXMLWriter->addInstance();
00889 }
00890
00891
00892 G4String hitType = "Hits";
00893 if (hitAttValues) {
00894 G4bool found = false;
00895 for (iAttVal = hitAttValues->begin();
00896 iAttVal != hitAttValues->end() && !found; ++iAttVal) {
00897 if (strcmp(iAttVal->GetName(),"HitType")==0) {
00898 hitType = iAttVal->GetValue();
00899 found = true;
00900 }
00901 }
00902 }
00903
00904
00905 G4String previousName = hepRepXMLWriter->prevTypeName[1];
00906 hepRepXMLWriter->addType(hitType,1);
00907
00908
00909
00910 if (strcmp(hitType,previousName)!=0) {
00911 hepRepXMLWriter->addAttValue("Layer",130);
00912
00913
00914
00915
00916
00917 if (hitAttValues && hitAttDefs) {
00918 for (iAttVal = hitAttValues->begin();
00919 iAttVal != hitAttValues->end(); ++iAttVal) {
00920 iAttDef = hitAttDefs->find(iAttVal->GetName());
00921 if (iAttDef != hitAttDefs->end()) {
00922
00923
00924 G4String category = iAttDef->second.GetCategory();
00925 if (strcmp(category,"Draw")!=0 &&
00926 strcmp(category,"Physics")!=0 &&
00927 strcmp(category,"Association")!=0 &&
00928 strcmp(category,"PickAction")!=0)
00929 category = "Physics";
00930 hepRepXMLWriter->addAttDef(iAttVal->GetName(), iAttDef->second.GetDesc(),
00931 category, iAttDef->second.GetExtra());
00932 }
00933 }
00934 }
00935 }
00936
00937
00938
00939 drawingHit = true;
00940 doneInitHit = false;
00941 G4VSceneHandler::AddCompound(hit);
00942 drawingHit = false;
00943 }
00944
00945
00946 void G4HepRepFileSceneHandler::InitTrajectory() {
00947 if (!doneInitTraj) {
00948
00949 hepRepXMLWriter->addInstance();
00950
00951
00952 if (trajAttValues) {
00953 std::vector<G4AttValue>::iterator iAttVal;
00954 for (iAttVal = trajAttValues->begin();
00955 iAttVal != trajAttValues->end(); ++iAttVal)
00956 hepRepXMLWriter->addAttValue(iAttVal->GetName(), iAttVal->GetValue());
00957 delete trajAttValues;
00958 }
00959
00960
00961 if (trajAttDefs)
00962 delete trajAttDefs;
00963
00964 doneInitTraj = true;
00965 }
00966 }
00967
00968
00969 void G4HepRepFileSceneHandler::InitHit() {
00970 if (!doneInitHit) {
00971
00972 hepRepXMLWriter->addInstance();
00973
00974
00975 if (hitAttValues) {
00976 std::vector<G4AttValue>::iterator iAttVal;
00977 for (iAttVal = hitAttValues->begin();
00978 iAttVal != hitAttValues->end(); ++iAttVal)
00979 hepRepXMLWriter->addAttValue(iAttVal->GetName(), iAttVal->GetValue());
00980 delete hitAttValues;
00981 }
00982
00983
00984 if (hitAttDefs)
00985 delete hitAttDefs;
00986
00987 doneInitHit = true;
00988 }
00989 }
00990
00991
00992 void G4HepRepFileSceneHandler::AddPrimitive(const G4Polyline& polyline) {
00993 #ifdef G4HEPREPFILEDEBUG
00994 G4cout <<
00995 "G4HepRepFileSceneHandler::AddPrimitive(const G4Polyline& polyline) called:"
00996 "\n polyline: " << polyline
00997 << G4endl;
00998 PrintThings();
00999 #endif
01000
01001 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01002
01003 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
01004 return;
01005
01006 if (inPrimitives2D) {
01007 if (!warnedAbout2DMarkers) {
01008 G4cout << "HepRepFile does not currently support 2D lines." << G4endl;
01009 warnedAbout2DMarkers = true;
01010 }
01011 return;
01012 }
01013
01014 if (drawingTraj)
01015 InitTrajectory();
01016
01017 if (drawingHit)
01018 InitHit();
01019
01020 haveVisible = true;
01021 AddHepRepInstance("Line", polyline);
01022
01023 hepRepXMLWriter->addPrimitive();
01024
01025 for (size_t i=0; i < polyline.size(); i++) {
01026 G4Point3D vertex = (fObjectTransformation) * polyline[i];
01027 hepRepXMLWriter->addPoint(vertex.x(), vertex.y(), vertex.z());
01028 }
01029 }
01030
01031
01032
01033 void G4HepRepFileSceneHandler::AddPrimitive (const G4Polymarker& line) {
01034 #ifdef G4HEPREPFILEDEBUG
01035 G4cout <<
01036 "G4HepRepFileSceneHandler::AddPrimitive(const G4Polymarker& line) called"
01037 << G4endl;
01038 PrintThings();
01039 #endif
01040
01041 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01042
01043 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
01044 return;
01045
01046 if (inPrimitives2D) {
01047 if (!warnedAbout2DMarkers) {
01048 G4cout << "HepRepFile does not currently support 2D lines." << G4endl;
01049 warnedAbout2DMarkers = true;
01050 }
01051 return;
01052 }
01053
01054 MarkerSizeType sizeType;
01055 G4double size = GetMarkerSize (line, sizeType);
01056 if (sizeType==world)
01057 size = 4.;
01058
01059 if (drawingTraj)
01060 return;
01061
01062 if (drawingHit)
01063 InitHit();
01064
01065 haveVisible = true;
01066 AddHepRepInstance("Point", line);
01067
01068 hepRepXMLWriter->addAttValue("MarkName", "Dot");
01069 hepRepXMLWriter->addAttValue("MarkSize", (G4int) size);
01070
01071 hepRepXMLWriter->addPrimitive();
01072
01073 for (size_t i=0; i < line.size(); i++) {
01074 G4Point3D vertex = (fObjectTransformation) * line[i];
01075 hepRepXMLWriter->addPoint(vertex.x(), vertex.y(), vertex.z());
01076 }
01077 }
01078
01079
01080 void G4HepRepFileSceneHandler::AddPrimitive(const G4Text& text) {
01081 #ifdef G4HEPREPFILEDEBUG
01082 G4cout <<
01083 "G4HepRepFileSceneHandler::AddPrimitive(const G4Text& text) called:"
01084 "\n text: " << text.GetText()
01085 << G4endl;
01086 PrintThings();
01087 #endif
01088
01089 if (!inPrimitives2D) {
01090 if (!warnedAbout3DText) {
01091 G4cout << "HepRepFile does not currently support 3D text." << G4endl;
01092 G4cout << "HepRep browsers can directly display text attributes on request." << G4endl;
01093 G4cout << "See Application Developers Guide for how to attach attributes to viewable objects." << G4endl;
01094 warnedAbout3DText = true;
01095 }
01096 return;
01097 }
01098
01099 MarkerSizeType sizeType;
01100 G4double size = GetMarkerSize (text, sizeType);
01101 if (sizeType==world)
01102 size = 12.;
01103
01104 haveVisible = true;
01105 AddHepRepInstance("Text", text);
01106
01107 hepRepXMLWriter->addAttValue("VAlign", "Top");
01108 hepRepXMLWriter->addAttValue("HAlign", "Left");
01109 hepRepXMLWriter->addAttValue("FontName", "Arial");
01110 hepRepXMLWriter->addAttValue("FontStyle", "Plain");
01111 hepRepXMLWriter->addAttValue("FontSize", (G4int) size);
01112 hepRepXMLWriter->addAttValue("FontHasBanner", "TRUE");
01113 hepRepXMLWriter->addAttValue("FontBannerColor", "0,0,0");
01114
01115 const G4Colour& colour = GetTextColour(text);
01116 float redness = colour.GetRed();
01117 float greenness = colour.GetGreen();
01118 float blueness = colour.GetBlue();
01119
01120
01121 if (redness==0. && greenness==0. && blueness==0.) {
01122 redness = 1.;
01123 greenness = 1.;
01124 blueness = 1.;
01125 }
01126 hepRepXMLWriter->addAttValue("FontColor",redness,greenness,blueness);
01127
01128 hepRepXMLWriter->addPrimitive();
01129
01130 hepRepXMLWriter->addAttValue("Text", text.GetText());
01131 hepRepXMLWriter->addAttValue("VPos", .99-text.GetYOffset());
01132 hepRepXMLWriter->addAttValue("HPos", text.GetXOffset());
01133 }
01134
01135
01136 void G4HepRepFileSceneHandler::AddPrimitive(const G4Circle& circle) {
01137 #ifdef G4HEPREPFILEDEBUG
01138 G4cout <<
01139 "G4HepRepFileSceneHandler::AddPrimitive(const G4Circle& circle) called:"
01140 "\n radius: " << circle.GetWorldRadius()
01141 << G4endl;
01142 PrintThings();
01143 #endif
01144
01145 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01146
01147 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
01148 return;
01149
01150 if (inPrimitives2D) {
01151 if (!warnedAbout2DMarkers) {
01152 G4cout << "HepRepFile does not currently support 2D circles." << G4endl;
01153 warnedAbout2DMarkers = true;
01154 }
01155 return;
01156 }
01157
01158 MarkerSizeType sizeType;
01159 G4double size = GetMarkerSize (circle, sizeType);
01160 if (sizeType==world)
01161 size = 4.;
01162
01163 if (drawingTraj)
01164 return;
01165
01166 if (drawingHit)
01167 InitHit();
01168
01169 haveVisible = true;
01170 AddHepRepInstance("Point", circle);
01171
01172 hepRepXMLWriter->addAttValue("MarkName", "Dot");
01173 hepRepXMLWriter->addAttValue("MarkSize", (G4int) size);
01174
01175 hepRepXMLWriter->addPrimitive();
01176
01177 G4Point3D center = (fObjectTransformation) * circle.GetPosition();
01178 hepRepXMLWriter->addPoint(center.x(), center.y(), center.z());
01179 }
01180
01181
01182 void G4HepRepFileSceneHandler::AddPrimitive(const G4Square& square) {
01183 #ifdef G4HEPREPFILEDEBUG
01184 G4cout <<
01185 "G4HepRepFileSceneHandler::AddPrimitive(const G4Square& square) called:"
01186 "\n side: " << square.GetWorldRadius()
01187 << G4endl;
01188 PrintThings();
01189 #endif
01190
01191 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01192
01193 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
01194 return;
01195
01196 if (inPrimitives2D) {
01197 if (!warnedAbout2DMarkers) {
01198 G4cout << "HepRepFile does not currently support 2D squares." << G4endl;
01199 warnedAbout2DMarkers = true;
01200 }
01201 return;
01202 }
01203
01204 MarkerSizeType sizeType;
01205 G4double size = GetMarkerSize (square, sizeType);
01206 if (sizeType==world)
01207 size = 4.;
01208
01209 if (drawingTraj)
01210 return;
01211
01212 if (drawingHit)
01213 InitHit();
01214
01215 haveVisible = true;
01216 AddHepRepInstance("Point", square);
01217
01218 hepRepXMLWriter->addAttValue("MarkName", "Square");
01219 hepRepXMLWriter->addAttValue("MarkSize", (G4int) size);
01220
01221 hepRepXMLWriter->addPrimitive();
01222
01223 G4Point3D center = (fObjectTransformation) * square.GetPosition();
01224 hepRepXMLWriter->addPoint(center.x(), center.y(), center.z());
01225 }
01226
01227
01228 void G4HepRepFileSceneHandler::AddPrimitive(const G4Polyhedron& polyhedron) {
01229 #ifdef G4HEPREPFILEDEBUG
01230 G4cout <<
01231 "G4HepRepFileSceneHandler::AddPrimitive(const G4Polyhedron& polyhedron) called."
01232 << G4endl;
01233 PrintThings();
01234 #endif
01235
01236 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01237
01238 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
01239 return;
01240
01241 if(polyhedron.GetNoFacets()==0)return;
01242
01243 if (drawingTraj)
01244 return;
01245
01246 if (drawingHit)
01247 InitHit();
01248
01249 haveVisible = true;
01250 AddHepRepInstance("Polygon", polyhedron);
01251
01252 G4Normal3D surfaceNormal;
01253 G4Point3D vertex;
01254
01255 G4bool notLastFace;
01256 do {
01257 hepRepXMLWriter->addPrimitive();
01258 notLastFace = polyhedron.GetNextNormal (surfaceNormal);
01259
01260 G4int edgeFlag = 1;
01261 G4bool notLastEdge;
01262 do {
01263 notLastEdge = polyhedron.GetNextVertex (vertex, edgeFlag);
01264 vertex = (fObjectTransformation) * vertex;
01265 hepRepXMLWriter->addPoint(vertex.x(), vertex.y(), vertex.z());
01266 } while (notLastEdge);
01267 } while (notLastFace);
01268 }
01269
01270
01271 void G4HepRepFileSceneHandler::AddPrimitive(const G4NURBS&) {
01272 #ifdef G4HEPREPFILEDEBUG
01273 G4cout <<
01274 "G4HepRepFileSceneHandler::AddPrimitive(const G4NURBS& nurbs) called."
01275 << G4endl;
01276 PrintThings();
01277 #endif
01278 G4cout << "G4HepRepFileSceneHandler::AddPrimitive G4NURBS : not implemented. " << G4endl;
01279 }
01280
01281
01282 G4HepRepFileXMLWriter *G4HepRepFileSceneHandler::GetHepRepXMLWriter() {
01283 return hepRepXMLWriter;
01284 }
01285
01286
01287 void G4HepRepFileSceneHandler::AddHepRepInstance(const char* primName,
01288 const G4Visible visible) {
01289 #ifdef G4HEPREPFILEDEBUG
01290 G4cout <<
01291 "G4HepRepFileSceneHandler::AddHepRepInstance called."
01292 << G4endl;
01293 #endif
01294
01295
01296 CheckFileOpen();
01297
01298 G4VPhysicalVolume* pCurrentPV = 0;
01299 G4LogicalVolume* pCurrentLV = 0;
01300 G4int currentDepth = 0;
01301 G4PhysicalVolumeModel* pPVModel = dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
01302 if (pPVModel) {
01303 pCurrentPV = pPVModel->GetCurrentPV();
01304 pCurrentLV = pPVModel->GetCurrentLV();
01305 currentDepth = pPVModel->GetCurrentDepth();
01306 }
01307
01308 #ifdef G4HEPREPFILEDEBUG
01309 G4cout <<
01310 "pCurrentPV:" << pCurrentPV << ", readyForTransients:" << fReadyForTransients
01311 << G4endl;
01312 #endif
01313
01314 if (drawingTraj || drawingHit) {
01315
01316
01317 }
01318 else if (fReadyForTransients) {
01319 if (strcmp("Event Data",hepRepXMLWriter->prevTypeName[0])!=0) {
01320 hepRepXMLWriter->addType("Event Data",0);
01321 hepRepXMLWriter->addInstance();
01322 }
01323
01324
01325
01326
01327
01328
01329
01330
01331
01332
01333
01334
01335
01336
01337
01338 int layer;
01339
01340 if (strcmp("Text",primName)==0) {
01341 hepRepXMLWriter->addType("EventID",1);
01342 } else {
01343 if (strcmp("Line",primName)==0) {
01344 hepRepXMLWriter->addType("TransientPolylines",1);
01345 layer = 100;
01346 } else {
01347 if (strcmp(hepRepXMLWriter->prevTypeName[1],"TransientPolylines")==0 &&
01348 strcmp("Square",primName)==0)
01349 {
01350 hepRepXMLWriter->addType("AuxiliaryPoints",2);
01351 layer = 110;
01352 } else {
01353 if (strcmp(hepRepXMLWriter->prevTypeName[1],"TransientPolylines")==0 &&
01354 strcmp("Circle",primName)==0)
01355 {
01356 hepRepXMLWriter->addType("StepPoints",2);
01357 layer = 120;
01358 } else {
01359 hepRepXMLWriter->addType("Hits",1);
01360 layer = 130;
01361 }
01362 }
01363 }
01364 hepRepXMLWriter->addAttValue("Layer",layer);
01365 }
01366
01367 hepRepXMLWriter->addInstance();
01368
01369
01370 } else if (pCurrentPV==0) {
01371 if (strcmp("AxesEtc",hepRepXMLWriter->prevTypeName[0])!=0) {
01372 hepRepXMLWriter->addType("AxesEtc",0);
01373 hepRepXMLWriter->addInstance();
01374 }
01375
01376 int layer;
01377
01378 if (strcmp("Text",primName)==0) {
01379 hepRepXMLWriter->addType("Text",1);
01380 } else {
01381 if (strcmp("Line",primName)==0) {
01382 hepRepXMLWriter->addType("Polylines",1);
01383 layer = 100;
01384 } else {
01385 hepRepXMLWriter->addType("Points",1);
01386 layer = 130;
01387 }
01388 hepRepXMLWriter->addAttValue("Layer",layer);
01389 }
01390
01391 hepRepXMLWriter->addInstance();
01392
01393
01394
01395
01396 } else {
01397
01398
01399 if (strcmp("Detector Geometry",hepRepXMLWriter->prevTypeName[0])!=0) {
01400
01401 hepRepXMLWriter->addType("Detector Geometry",0);
01402 hepRepXMLWriter->addInstance();
01403 }
01404
01405
01406
01407 if(strcmp(hepRepXMLWriter->prevTypeName[currentDepth+1],pCurrentPV->GetName())!=0) {
01408
01409 typedef G4PhysicalVolumeModel::G4PhysicalVolumeNodeID PVNodeID;
01410 typedef std::vector<PVNodeID> PVPath;
01411 const PVPath& drawnPVPath = pPVModel->GetDrawnPVPath();
01412 PVPath::const_reverse_iterator ri = ++drawnPVPath.rbegin();
01413 G4int drawnMotherDepth;
01414 if (ri != drawnPVPath.rend()) {
01415
01416 drawnMotherDepth = ri->GetNonCulledDepth();
01417
01418 } else {
01419
01420 drawnMotherDepth = -1;
01421
01422 }
01423
01424 while (drawnMotherDepth < (currentDepth-1)) {
01425 G4String culledParentName = "Culled parent of " + pCurrentPV->GetName();
01426
01427 hepRepXMLWriter->addType(culledParentName, drawnMotherDepth+2);
01428 hepRepXMLWriter->addInstance();
01429 drawnMotherDepth ++;
01430 }
01431 }
01432
01433
01434 hepRepXMLWriter->addType(pCurrentPV->GetName(),currentDepth+1);
01435 hepRepXMLWriter->addInstance();
01436
01437 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01438
01439 if (fpVisAttribs && (fpVisAttribs->IsVisible()==0) && messenger->getCullInvisibles())
01440 return;
01441
01442
01443 hepRepXMLWriter->addAttValue("Layer",hepRepXMLWriter->typeDepth);
01444 hepRepXMLWriter->addAttValue("LVol", pCurrentLV->GetName());
01445 G4Region* region = pCurrentLV->GetRegion();
01446 G4String regionName = region? region->GetName(): G4String("No region");
01447 hepRepXMLWriter->addAttValue("Region", regionName);
01448 hepRepXMLWriter->addAttValue("RootRegion", pCurrentLV->IsRootRegion());
01449 hepRepXMLWriter->addAttValue("Solid", pCurrentLV->GetSolid()->GetName());
01450 hepRepXMLWriter->addAttValue("EType", pCurrentLV->GetSolid()->GetEntityType());
01451 G4Material * material = pPVModel->GetCurrentMaterial();
01452 G4String matName = material? material->GetName(): G4String("No material");
01453 hepRepXMLWriter->addAttValue("Material", matName);
01454 G4double matDensity = material? material->GetDensity(): 0.;
01455 hepRepXMLWriter->addAttValue("Density", matDensity*m3/kg);
01456 G4State matState = material? material->GetState(): kStateUndefined;
01457 hepRepXMLWriter->addAttValue("State", matState);
01458 G4double matRadlen = material? material->GetRadlen(): 0.;
01459 hepRepXMLWriter->addAttValue("Radlen", matRadlen/m);
01460 }
01461
01462 hepRepXMLWriter->addAttValue("DrawAs",primName);
01463
01464
01465 float redness;
01466 float greenness;
01467 float blueness;
01468 G4bool isVisible;
01469
01470 if (fpVisAttribs || haveVisible) {
01471 G4Colour colour;
01472
01473 if (fpVisAttribs) {
01474 colour = fpVisAttribs->GetColour();
01475 isVisible = fpVisAttribs->IsVisible();
01476 } else {
01477 colour = GetColour(visible);
01478 isVisible = fpViewer->
01479 GetApplicableVisAttributes(visible.GetVisAttributes())->IsVisible();
01480 }
01481
01482 redness = colour.GetRed();
01483 greenness = colour.GetGreen();
01484 blueness = colour.GetBlue();
01485
01486
01487 if (redness==0. && greenness==0. && blueness==0.) {
01488 redness = 1.;
01489 greenness = 1.;
01490 blueness = 1.;
01491 }
01492 } else {
01493 #ifdef G4HEPREPFILEDEBUG
01494 G4cout <<
01495 "G4HepRepFileSceneHandler::AddHepRepInstance using default colour."
01496 << G4endl;
01497 #endif
01498 redness = 1.;
01499 greenness = 1.;
01500 blueness = 1.;
01501 isVisible = true;
01502 }
01503
01504 if (strcmp(primName,"Point")==0)
01505 hepRepXMLWriter->addAttValue("MarkColor",redness,greenness,blueness);
01506 else
01507 hepRepXMLWriter->addAttValue("LineColor",redness,greenness,blueness);
01508
01509 hepRepXMLWriter->addAttValue("Visibility",isVisible);
01510 }
01511
01512
01513 void G4HepRepFileSceneHandler::CheckFileOpen() {
01514 #ifdef G4HEPREPFILEDEBUG
01515 G4cout <<
01516 "G4HepRepFileSceneHandler::CheckFileOpen called."
01517 << G4endl;
01518 #endif
01519
01520 if (!hepRepXMLWriter->isOpen) {
01521 G4String newFileSpec;
01522
01523 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01524
01525 if (messenger->getOverwrite()) {
01526 newFileSpec = messenger->getFileDir()+messenger->getFileName()+".heprep";
01527 } else {
01528 newFileSpec = messenger->getFileDir()+messenger->getFileName()+G4UIcommand::ConvertToString(fileCounter)+".heprep";
01529 }
01530
01531 G4cout << "HepRepFile writing to " << newFileSpec << G4endl;
01532
01533 hepRepXMLWriter->open(newFileSpec);
01534
01535 if (!messenger->getOverwrite())
01536 fileCounter++;
01537
01538 hepRepXMLWriter->addAttDef("Generator", "HepRep Data Generator", "Physics","");
01539 G4String versionString = G4Version;
01540 versionString = versionString.substr(1,versionString.size()-2);
01541 versionString = " Geant4 version " + versionString + " " + G4Date;
01542 hepRepXMLWriter->addAttValue("Generator", versionString);
01543
01544 hepRepXMLWriter->addAttDef("LVol", "Logical Volume", "Physics","");
01545 hepRepXMLWriter->addAttDef("Region", "Cuts Region", "Physics","");
01546 hepRepXMLWriter->addAttDef("RootRegion", "Root Region", "Physics","");
01547 hepRepXMLWriter->addAttDef("Solid", "Solid Name", "Physics","");
01548 hepRepXMLWriter->addAttDef("EType", "Entity Type", "Physics","");
01549 hepRepXMLWriter->addAttDef("Material", "Material Name", "Physics","");
01550 hepRepXMLWriter->addAttDef("Density", "Material Density", "Physics","kg/m3");
01551 hepRepXMLWriter->addAttDef("State", "Material State", "Physics","");
01552 hepRepXMLWriter->addAttDef("Radlen", "Material Radiation Length", "Physics","m");
01553 }
01554 }
01555
01556
01557 void G4HepRepFileSceneHandler::ClearTransientStore() {
01558
01559
01560
01561 if (fpViewer) {
01562 fpViewer -> SetView();
01563 fpViewer -> ClearView();
01564 fpViewer -> DrawView();
01565 }
01566 }