G4HepRepSceneHandler.cc

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration.                      *
00020 // * By using,  copying,  modifying or  distributing the software (or *
00021 // * any work based  on the software)  you  agree  to acknowledge its *
00022 // * use  in  resulting  scientific  publications,  and indicate your *
00023 // * acceptance of all terms of the Geant4 Software license.          *
00024 // ********************************************************************
00025 //
00026 // $Id$
00027 //
00028 
00033 #include <stdio.h>
00034 
00035 #include "globals.hh"
00036 #include <vector>
00037 #include <iostream>
00038 // NOTE not available on Solaris 5.2 and Linux g++ 2.95.2
00039 // #include <sstream>
00040 #include <iomanip>
00041 #include <fstream>
00042 #include <cmath>
00043 
00044 //HepRep
00045 #include "HEPREP/HepRep.h"
00046 #include "G4HepRepMessenger.hh"
00047 
00048 //G4
00049 #include "G4PhysicalConstants.hh"
00050 #include "G4Vector3D.hh"
00051 #include "G4Version.hh"
00052 #include "G4Types.hh"
00053 #include "G4Point3D.hh"
00054 #include "G4Normal3D.hh"
00055 #include "G4Polyline.hh"
00056 #include "G4Polymarker.hh"
00057 #include "G4Polyhedron.hh"
00058 #include "G4Circle.hh"
00059 #include "G4Square.hh"
00060 #include "G4Text.hh"
00061 #include "G4NURBS.hh"
00062 #include "G4VPhysicalVolume.hh"
00063 #include "G4VisAttributes.hh"
00064 #include "G4VSolid.hh"
00065 #include "G4VTrajectory.hh"
00066 #include "G4VTrajectoryPoint.hh"
00067 #include "G4VHit.hh"
00068 #include "G4Scene.hh"
00069 #include "G4Material.hh"
00070 #include "G4AttDef.hh"
00071 #include "G4AttValue.hh"
00072 #include "G4AttCheck.hh"
00073 
00074 // CHepRep
00075 #include "cheprep/XMLHepRepFactory.h"
00076 
00077 // This
00078 #include "G4HepRep.hh"
00079 #include "G4HepRepSceneHandler.hh"
00080 #include "G4HepRepViewer.hh"
00081 
00082 
00083 using namespace HEPREP;
00084 using namespace cheprep;
00085 using namespace std;
00086 
00087 G4int G4HepRepSceneHandler::sceneIdCount = 0;
00088 
00089 //#define LDEBUG 1
00090 //#define SDEBUG 1
00091 //#define PDEBUG 1
00092 
00093 G4HepRepSceneHandler::G4HepRepSceneHandler (G4VGraphicsSystem& system, const G4String& name)
00094         : G4VSceneHandler (system, sceneIdCount++, name),
00095           geometryLayer         ("Geometry"),
00096           eventLayer            ("Event"),
00097           calHitLayer           ("CalHit"),
00098           trajectoryLayer       ("Trajectory"),
00099           hitLayer              ("Hit"),
00100           rootVolumeName        ("Geometry"),
00101           baseName              (""),
00102           eventNumberPrefix     (""),
00103           eventNumberSuffix     (""),
00104           eventNumber           (1),
00105           eventNumberWidth      (-1),
00106           extension             (""),
00107           writeBinary           (false),
00108           writeZip              (false),
00109           writeGZ               (false),
00110           writeMultipleFiles    (false),
00111           currentHit            (NULL),
00112           currentTrack          (NULL),
00113           _heprep               (NULL),
00114           _heprepGeometry       (NULL)
00115 {
00116 
00117 #ifdef LDEBUG
00118     cout << "G4HepRepSceneHandler::G4HepRepSceneHandler: " << system << endl;
00119 #endif
00120 
00121     materialState[kStateSolid]      = G4String("Solid");
00122     materialState[kStateLiquid]     = G4String("Liquid");
00123     materialState[kStateGas]        = G4String("Gas");
00124     materialState[kStateUndefined]  = G4String("Undefined");    
00125 
00126     factory = new XMLHepRepFactory();
00127     writer = NULL;
00128     
00129     // opening of file deferred to closeHepRep();
00130     openHepRep();
00131 }
00132 
00133 
00134 G4HepRepSceneHandler::~G4HepRepSceneHandler () {
00135 #ifdef LDEBUG
00136     cout << "G4HepRepSceneHandler::~G4HepRepSceneHandler() " << endl;
00137 #endif
00138     close();
00139 
00140     delete factory;
00141     factory = NULL;
00142 
00143     dynamic_cast<G4HepRep*>(GetGraphicsSystem())->removeSceneHandler();
00144 }
00145 
00146 
00147 void G4HepRepSceneHandler::open(G4String name) {
00148     if (writer != NULL) return;
00149 
00150     if (name == "stdout") {
00151 #ifdef LDEBUG
00152         cout << "G4HepRepSceneHandler::Open() stdout" << endl;
00153 #endif
00154         writer = factory->createHepRepWriter(&cout, false, false);
00155         out  = NULL;
00156         baseName = name;
00157         eventNumberPrefix = "";
00158         eventNumberSuffix = "";
00159         extension = "";
00160         writeBinary = false;
00161         writeZip = false;
00162         writeGZ = false;
00163         writeMultipleFiles = false;
00164         eventNumber = 0;
00165         eventNumberWidth = 0;        
00166     } else if (name == "stderr") {
00167 #ifdef LDEBUG
00168         cout << "G4HepRepSceneHandler::Open() stderr" << endl;
00169 #endif
00170         writer = factory->createHepRepWriter(&cerr, false, false);
00171         out = NULL;
00172         baseName = name;
00173         eventNumberPrefix = "";
00174         eventNumberSuffix = "";
00175         extension = "";
00176         writeBinary = false;
00177         writeZip = false;
00178         writeGZ = false;
00179         writeMultipleFiles = false;
00180         eventNumber = 0;
00181         eventNumberWidth = 0;        
00182     } else {
00183 #ifdef LDEBUG
00184         cout << "G4HepRepSceneHandler::Open() " << name << endl;
00185 #endif
00186         if (eventNumberWidth < 0) {
00187             // derive filename(s)
00188             // check for extensions
00189             const unsigned int numberOfExtensions = 8;
00190             string ext[numberOfExtensions] = {".heprep", ".heprep.xml", ".heprep.zip", ".heprep.gz",
00191                                               ".bheprep", ".bheprep.xml", ".bheprep.zip", ".bheprep.gz"};
00192             unsigned int i=0;
00193             while (i < numberOfExtensions) {
00194                 int dot = name.size() - ext[i].size();
00195                 if ((dot >= 0) && 
00196                     (name.substr(dot, ext[i].size()) == ext[i])) break;
00197                 i++;
00198             }
00199         
00200             if (i != numberOfExtensions) {
00201                 extension = ext[i];
00202                 writeBinary = i >= (numberOfExtensions/2);
00203                 writeZip = (i == 2) || (i == 6);
00204                 writeGZ = (i == 3) || (i == 7);
00205 
00206                 int dot = name.length() - extension.length();
00207                 baseName = (dot >= 0) ? name.substr(0, dot) : "";
00208 
00209 #ifndef G4LIB_USE_ZLIB               
00210                 if (writeGZ) {
00211                     cerr << endl;
00212                     cerr << "WARNING: the .gz output file you are creating will be a plain file," << endl;
00213                     cerr << "       since compression support (ZLIB) was not compiled into the Geant4." << endl;
00214                     cerr << "       To avoid confusion with real gz files, the output filename has been" << endl;
00215                     cerr << "       extended with the name '.no-gz'." << endl;
00216                     cerr << "       A plain heprep or bheprep file can be fairly large." << endl; 
00217                 }
00218                 if (writeZip) {
00219                     cerr << endl;
00220                     cerr << "WARNING: the .zip output file you are creating will not be compressed," << endl;
00221                     cerr << "       since compression support (ZLIB) was not compiled into the Geant4." << endl;
00222                     cerr << "       A zip file containing non-compressed heprep or bheprep files can" << endl;
00223                     cerr << "       be fairly large." << endl;
00224                 }
00225                 if (writeGZ || writeZip) {
00226                     cerr << "SOLUTION: To add compression support using ZLIB, you need to:" << endl;
00227                     cerr << "       1. Define G4LIB_USE_ZLIB and recompile the visualization category." << endl;
00228                     cerr << "       2. Optionally define G4_LIB_BUILD_ZLIB if your system does not have" << endl;
00229                     cerr << "          zlib installed (e.g. WIN32-VC)." << endl;
00230                     cerr << "       3. Relink your application code." << endl;  
00231                     cerr << endl;
00232                 }
00233                 if (writeGZ) {
00234                     extension = extension + ".no-gz";
00235                     writeGZ = false;
00236                 }
00237 #endif // G4LIB_USE_ZLIB
00238             } else {
00239                 // Default for no extension  
00240                 extension = ".heprep.zip";
00241                 writeBinary = false;
00242                 writeZip = true;
00243                 writeGZ = false;
00244                 baseName = name;
00245             }
00246         
00247             writeMultipleFiles = false;
00248             int startDigit = -1; int endDigit = -1;
00249 
00250                         G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00251 
00252             string suffix = messenger->getEventNumberSuffix();
00253             if (suffix != "") {
00254                 // look for 0000 pattern in suffix
00255                 endDigit = suffix.length()-1; 
00256                 while (endDigit >= 0) {
00257                     if (isdigit(suffix.at(endDigit))) break;
00258                     endDigit--;
00259                 }
00260                 if (endDigit < 0) {
00261                     cerr << "/vis/heprep/appendEventNumberSuffix contains no digits" << endl;
00262                 } else {
00263                     writeMultipleFiles = true;
00264                     startDigit = endDigit;
00265                     while (startDigit >= 0) {
00266                         if (!isdigit(suffix.at(startDigit))) break;
00267                         startDigit--;
00268                     }
00269                     startDigit++;
00270                 }                
00271             }
00272 
00273             if (writeMultipleFiles) {
00274                 eventNumberPrefix = suffix.substr(0, startDigit);
00275                 eventNumber = atoi(suffix.substr(startDigit, endDigit).c_str());
00276                 eventNumberWidth = endDigit +1 - startDigit;
00277                 eventNumberSuffix = suffix.substr(endDigit+1);
00278             } else {
00279                 // open single file here
00280                 openFile(baseName+extension);
00281         
00282                 eventNumber = 1;
00283                 eventNumberWidth = 10;
00284                 eventNumberPrefix = "";
00285                 eventNumberSuffix = "";
00286             }   
00287         }    
00288     }
00289 }
00290 
00291 
00292 void G4HepRepSceneHandler::openHepRep() {
00293 #ifdef LDEBUG
00294     cout << "G4HepRepSceneHandler::OpenHepRep() " << endl;
00295 #endif
00296 
00297     if (_heprep != NULL) return;
00298 
00299     // all done on demand, once pointers are set to NULL
00300     _heprepGeometry       = NULL;
00301     _geometryInstanceTree = NULL;
00302     _geometryRootInstance = NULL;
00303     _geometryInstance.clear();
00304     _geometryTypeTree     = NULL;
00305     _geometryRootType     = NULL;
00306     _geometryTypeName.clear();
00307     _geometryType.clear();
00308     _eventInstanceTree    = NULL;
00309     _eventInstance        = NULL;
00310     _eventTypeTree        = NULL;
00311     _eventType            = NULL;
00312     _trajectoryType       = NULL;
00313     _hitType              = NULL;
00314     _calHitType           = NULL;
00315     _calHitFaceType       = NULL;     
00316 }
00317 
00318 
00322 bool G4HepRepSceneHandler::closeHepRep(bool final) {
00323     if (_heprep == NULL) return true;
00324 
00325 #ifdef LDEBUG
00326     cout << "G4HepRepSceneHandler::CloseHepRep() start" << endl;
00327 #endif
00328 
00329     // if this is the final close, then there should not be any event pending to be written.
00330     if (final) {
00331         if (_eventInstanceTree != NULL) {
00332             cerr << "WARNING: you probably used '/vis/viewer/endOfEventAction accumulate' and "
00333                  << "forgot to call /vis/viewer/update before exit. No event written." << endl;
00334         }
00335     } else {
00336                 
00337                 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00338 
00339         // add geometry to the heprep if there is an event (separate geometries are written
00340         // using DrawView() called from /vis/viewer/flush)
00341         if (_eventInstanceTree != NULL) {
00342             GetCurrentViewer()->DrawView();
00343         
00344             // couple geometry
00345 
00346             if ( messenger->appendGeometry()) {
00347                 // couple geometry to event if geometry was written
00348                 if ((_geometryInstanceTree != NULL)) {
00349                     getEventInstanceTree()->addInstanceTree(getGeometryInstanceTree());
00350                 }
00351             } else {
00352                 char name[128];
00353                 if (writeMultipleFiles) {
00354                     sprintf(name, "%s%s%s#%s", baseName.c_str(), "-geometry", extension.c_str(), "G4GeometryData");
00355                 } else {
00356                     sprintf(name, "%s%s#%s", "geometry", (writeBinary ? ".bheprep" : ".heprep"), "G4GeometryData");
00357                 }   
00358                 getEventInstanceTree()->addInstanceTree(factory->createHepRepTreeID(name, "1.0"));
00359             }
00360         }
00361     
00362         // force inclusion of all subtypes of event
00363         if (_eventInstanceTree != NULL) {
00364             getEventType();
00365             getTrajectoryType();
00366             getHitType();
00367             getCalHitType();
00368             getCalHitFaceType();
00369         }
00370     
00371         // Give this HepRep all of the layer order info for both geometry and event,
00372         // since these will both end up in a single HepRep.
00373         writeLayers(_heprepGeometry);
00374         writeLayers(_heprep);
00375 
00376         // open heprep file
00377         if (writer == NULL) {
00378             open((GetScene() == NULL) ? G4String("G4HepRepOutput.heprep.zip") : GetScene()->GetName());
00379         }
00380 
00381         // write out separate geometry
00382         if (! messenger->appendGeometry() && (_heprepGeometry != NULL)) {
00383             if (writeMultipleFiles) {
00384                 char fileName[128];
00385                 sprintf(fileName, "%s%s%s", baseName.c_str(), "-geometry", extension.c_str());
00386                 openFile(G4String(fileName));
00387             }
00388 
00389             char name[128];
00390             sprintf(name, "%s%s", "geometry", (writeBinary ? ".bheprep" : ".heprep"));
00391             if (!writeMultipleFiles) {
00392                 writer->addProperty("RecordLoop.ignore", name);
00393             }
00394 
00395             writer->write(_heprepGeometry, G4String(name));
00396             
00397             delete _heprepGeometry;
00398             _heprepGeometry = NULL;
00399             
00400             if (writeMultipleFiles) closeFile();
00401         }
00402     
00403         if (writeMultipleFiles) {
00404 // NOTE: does not work on Solaris 5.2 and Linux 2.95.2
00405 //            stringstream fileName;
00406 //            fileName << baseName << eventNumberPrefix << setw(eventNumberWidth) << setfill('0') << eventNumber << eventNumberSuffix << extension;
00407 //            openFile(fileName.str());
00408 // Use instead:
00409             char fileName[128];
00410             char fileFormat[128];
00411             sprintf(fileFormat, "%s%d%s", "%s%s%0", eventNumberWidth, "d%s%s");
00412             sprintf(fileName, fileFormat, baseName.c_str(), eventNumberPrefix.c_str(), eventNumber, eventNumberSuffix.c_str(), extension.c_str());
00413             openFile(G4String(fileName));
00414         }
00415             
00416         // write out the heprep
00417 // NOTE: does not work on Solaris 5.2 and Linux 2.95.2
00418 //        stringstream eventName;
00419 //        eventName << "event-" << setw(eventNumberWidth) << setfill('0') << eventNumber << (writeBinary ? ".bheprep" : ".heprep");
00420 //        writer->write(_heprep, eventName.str());
00421 // Use instead:
00422         char eventName[128];
00423         char eventFormat[128];
00424         sprintf(eventFormat, "%s%d%s%s", "event-%0", eventNumberWidth, "d", (writeBinary ? ".bheprep" : ".heprep"));
00425         sprintf(eventName, eventFormat, eventNumber);
00426         writer->write(_heprep, G4String(eventName));
00427 
00428         eventNumber++;
00429     }
00430     
00431     delete _heprep;
00432     _heprep = NULL;
00433 
00434     if (writeMultipleFiles) closeFile();
00435 
00436     return true;
00437 }
00438 
00439 
00440 void G4HepRepSceneHandler::close() {
00441 
00442 #ifdef LDEBUG
00443     cout << "G4HepRepSceneHandler::Close() " << endl;
00444 #endif
00445 
00446     if (writer == NULL) return;
00447     
00448     if (!writeMultipleFiles) {
00449         closeHepRep(true);
00450         closeFile();
00451     }
00452     
00453     G4HepRepViewer* viewer = dynamic_cast<G4HepRepViewer*>(GetCurrentViewer());
00454     viewer->reset();
00455 }
00456 
00457 void G4HepRepSceneHandler::openFile(G4String name) {
00458     out = new ofstream(name.c_str(), std::ios::out | std::ios::binary );
00459     writer = factory->createHepRepWriter(out, writeZip, writeZip || writeGZ);
00460 }
00461 
00462 void G4HepRepSceneHandler::closeFile() {
00463     writer->close();
00464     delete writer;
00465     writer = NULL;
00466 
00467     delete out;
00468     out = NULL;
00469 }
00470 
00471 void G4HepRepSceneHandler::writeLayers(HepRep* heprep) {
00472     if (heprep == NULL) return;
00473     heprep->addLayer(geometryLayer); 
00474     heprep->addLayer(eventLayer);
00475     heprep->addLayer(calHitLayer);
00476     heprep->addLayer(trajectoryLayer);
00477     heprep->addLayer(hitLayer);
00478 }  
00479 
00480 void G4HepRepSceneHandler::BeginModeling() {
00481 #ifdef SDEBUG
00482     cout << "G4HepRepSceneHandler::BeginModeling() " << endl;
00483 #endif
00484     G4VSceneHandler::BeginModeling();
00485 }
00486 
00487 
00488 void G4HepRepSceneHandler::EndModeling() {
00489 #ifdef SDEBUG
00490     cout << "G4HepRepSceneHandler::EndModeling() " << endl;
00491 #endif
00492     G4VSceneHandler::EndModeling();
00493 }
00494 
00495 void G4HepRepSceneHandler::AddSolid(const G4Box& box) {
00496 #ifdef SDEBUG
00497     cout << "G4HepRepSceneHandler::AddSolid(const G4Box& box)" << endl;
00498 #endif
00499 
00500     if (dontWrite()) return;
00501         
00502         G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00503 
00504     if (! messenger->useSolids()) {
00505         G4VSceneHandler::AddSolid(box);
00506         return;
00507     }
00508     
00509     G4double dx = box.GetXHalfLength();
00510     G4double dy = box.GetYHalfLength();
00511     G4double dz = box.GetZHalfLength();
00512     
00513     G4Point3D vertex1(G4Point3D( dx, dy,-dz));
00514     G4Point3D vertex2(G4Point3D( dx,-dy,-dz));
00515     G4Point3D vertex3(G4Point3D(-dx,-dy,-dz));
00516     G4Point3D vertex4(G4Point3D(-dx, dy,-dz));
00517     G4Point3D vertex5(G4Point3D( dx, dy, dz));
00518     G4Point3D vertex6(G4Point3D( dx,-dy, dz));
00519     G4Point3D vertex7(G4Point3D(-dx,-dy, dz));
00520     G4Point3D vertex8(G4Point3D(-dx, dy, dz));
00521     
00522     vertex1 = (transform) * vertex1;
00523     vertex2 = (transform) * vertex2;
00524     vertex3 = (transform) * vertex3;
00525     vertex4 = (transform) * vertex4;
00526     vertex5 = (transform) * vertex5;
00527     vertex6 = (transform) * vertex6;
00528     vertex7 = (transform) * vertex7;
00529     vertex8 = (transform) * vertex8;
00530 
00531     HepRepInstance* instance = getGeometryOrEventInstance(getCalHitType());
00532     addAttributes(instance, getCalHitType());
00533 
00534     setAttribute(instance, "DrawAs", G4String("Prism"));
00535         
00536     setVisibility(instance, box);
00537     setLine(instance, box);
00538     setColor(instance, getColorFor(box));
00539 
00540     factory->createHepRepPoint(instance, vertex1.x(), vertex1.y(), vertex1.z());
00541     factory->createHepRepPoint(instance, vertex2.x(), vertex2.y(), vertex2.z());
00542     factory->createHepRepPoint(instance, vertex3.x(), vertex3.y(), vertex3.z());
00543     factory->createHepRepPoint(instance, vertex4.x(), vertex4.y(), vertex4.z());
00544     factory->createHepRepPoint(instance, vertex5.x(), vertex5.y(), vertex5.z());
00545     factory->createHepRepPoint(instance, vertex6.x(), vertex6.y(), vertex6.z());
00546     factory->createHepRepPoint(instance, vertex7.x(), vertex7.y(), vertex7.z());
00547     factory->createHepRepPoint(instance, vertex8.x(), vertex8.y(), vertex8.z());
00548 }
00549 
00550 
00551 void G4HepRepSceneHandler::AddSolid(const G4Cons& cons) {
00552 #ifdef SDEBUG
00553     cout << "G4HepRepSceneHandler::AddSolid(const G4Cons& cons)" << endl;
00554 #endif
00555 
00556     if (dontWrite()) return;
00557         
00558         G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00559 
00560     if (! messenger->useSolids() || (cons.GetDeltaPhiAngle() < twopi)) {
00561         G4VSceneHandler::AddSolid(cons);
00562         return;
00563     }
00564 
00565     G4PhysicalVolumeModel* pPVModel =
00566       dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
00567     if (!pPVModel) {
00568       G4VSceneHandler::AddSolid(cons);
00569         return;
00570     }
00571 
00572     G4LogicalVolume* pCurrentLV = pPVModel->GetCurrentLV();
00573     G4int currentDepth = pPVModel->GetCurrentDepth();
00574     G4Material* pCurrentMaterial = pPVModel->GetCurrentMaterial();
00575 
00576     G4Point3D vertex1(G4Point3D( 0., 0., cons.GetZHalfLength()));
00577     G4Point3D vertex2(G4Point3D( 0., 0.,-cons.GetZHalfLength()));
00578 
00579     vertex1 = (transform) * vertex1;
00580     vertex2 = (transform) * vertex2;
00581 
00582     HepRepInstance* instance = getGeometryInstance(pCurrentLV, pCurrentMaterial, currentDepth);
00583     setAttribute(instance, "DrawAs", G4String("Cylinder"));
00584         
00585     setVisibility(instance, cons);
00586     setLine(instance, cons);
00587     setColor(instance, getColorFor(cons));
00588 
00589     HepRepType* type = getGeometryType(pCurrentLV->GetName(), currentDepth);
00590 
00591     // Outer cylinder.
00592     HepRepInstance* outer = factory->createHepRepInstance(instance, type);
00593     outer->addAttValue("pickParent",true);
00594     outer->addAttValue("showParentAttributes",true);
00595     
00596     HepRepPoint* op1 = factory->createHepRepPoint(outer, vertex1.x(), vertex1.y(), vertex1.z());
00597     op1->addAttValue("Radius",cons.GetOuterRadiusPlusZ());
00598     
00599     HepRepPoint* op2 = factory->createHepRepPoint(outer, vertex2.x(), vertex2.y(), vertex2.z());
00600     op2->addAttValue("Radius",cons.GetOuterRadiusMinusZ());
00601 
00602     // Inner cylinder.
00603     HepRepInstance* inner = factory->createHepRepInstance(instance, type);
00604     inner->addAttValue("pickParent",true);
00605     inner->addAttValue("showParentAttributes",true);
00606     
00607     HepRepPoint* ip1 = factory->createHepRepPoint(inner, vertex1.x(), vertex1.y(), vertex1.z());
00608     ip1->addAttValue("Radius",cons.GetInnerRadiusPlusZ());
00609     
00610     HepRepPoint* ip2 = factory->createHepRepPoint(inner, vertex2.x(), vertex2.y(), vertex2.z());
00611     ip2->addAttValue("Radius",cons.GetInnerRadiusMinusZ());
00612 }
00613 
00614 
00615 void G4HepRepSceneHandler::AddSolid(const G4Tubs& tubs) {
00616 #ifdef SDEBUG
00617     cout << "G4HepRepSceneHandler::AddSolid(const G4Tubs& tubs)" << endl;
00618 #endif
00619 
00620     if (dontWrite()) return;
00621         
00622         G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00623 
00624     if (! messenger->useSolids() || (tubs.GetDeltaPhiAngle() < twopi)) {
00625         G4VSceneHandler::AddSolid(tubs);
00626         return;
00627     }
00628     
00629     G4PhysicalVolumeModel* pPVModel =
00630       dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
00631     if (!pPVModel) {
00632       G4VSceneHandler::AddSolid(tubs);
00633         return;
00634     }
00635 
00636     G4LogicalVolume* pCurrentLV = pPVModel->GetCurrentLV();
00637     G4int currentDepth = pPVModel->GetCurrentDepth();
00638     G4Material* pCurrentMaterial = pPVModel->GetCurrentMaterial();
00639 
00640     G4Point3D vertex1(G4Point3D( 0., 0., tubs.GetZHalfLength()));
00641     G4Point3D vertex2(G4Point3D( 0., 0.,-tubs.GetZHalfLength()));
00642 
00643     vertex1 = (transform) * vertex1;
00644     vertex2 = (transform) * vertex2;
00645 
00646     HepRepInstance* instance = getGeometryInstance(pCurrentLV, pCurrentMaterial, currentDepth);
00647     setAttribute(instance, "DrawAs", G4String("Cylinder"));
00648         
00649     setVisibility(instance, tubs);
00650     setLine(instance, tubs);
00651     setColor(instance, getColorFor(tubs));
00652 
00653     HepRepType* type = getGeometryType(pCurrentLV->GetName(), currentDepth);
00654 
00655     // Outer cylinder.
00656     HepRepInstance* outer = factory->createHepRepInstance(instance, type);
00657     outer->addAttValue("Radius",tubs.GetOuterRadius());
00658     outer->addAttValue("pickParent",true);
00659     outer->addAttValue("showParentAttributes",true);
00660     factory->createHepRepPoint(outer, vertex1.x(), vertex1.y(), vertex1.z());
00661     factory->createHepRepPoint(outer, vertex2.x(), vertex2.y(), vertex2.z());
00662 
00663     // Inner cylinder.
00664     if (tubs.GetInnerRadius() > 0.) {
00665         HepRepInstance* inner = factory->createHepRepInstance(instance, type);
00666         inner->addAttValue("Radius",tubs.GetInnerRadius());
00667         inner->addAttValue("pickParent",true);
00668         inner->addAttValue("showParentAttributes",true);
00669         factory->createHepRepPoint(inner, vertex1.x(), vertex1.y(), vertex1.z());
00670         factory->createHepRepPoint(inner, vertex2.x(), vertex2.y(), vertex2.z());
00671     }
00672 }
00673 
00674 
00675 void G4HepRepSceneHandler::AddSolid(const G4Trd& trd) {
00676 #ifdef SDEBUG
00677     cout << "G4HepRepSceneHandler::AddSolid(const G4Trd& trd)" << endl;
00678 #endif
00679     if (dontWrite()) return;
00680         
00681         G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
00682 
00683     if (! messenger->useSolids()) {
00684         G4VSceneHandler::AddSolid(trd);
00685         return;
00686     }
00687     
00688     G4double dx1 = trd.GetXHalfLength1();
00689     G4double dy1 = trd.GetYHalfLength1();
00690     G4double dx2 = trd.GetXHalfLength2();
00691     G4double dy2 = trd.GetYHalfLength2();
00692     G4double dz = trd.GetZHalfLength();
00693     
00694     G4Point3D vertex1(G4Point3D( dx1, dy1,-dz));
00695     G4Point3D vertex2(G4Point3D( dx1,-dy1,-dz));
00696     G4Point3D vertex3(G4Point3D(-dx1,-dy1,-dz));
00697     G4Point3D vertex4(G4Point3D(-dx1, dy1,-dz));
00698     G4Point3D vertex5(G4Point3D( dx2, dy2, dz));
00699     G4Point3D vertex6(G4Point3D( dx2,-dy2, dz));
00700     G4Point3D vertex7(G4Point3D(-dx2,-dy2, dz));
00701     G4Point3D vertex8(G4Point3D(-dx2, dy2, dz));
00702     
00703     vertex1 = (transform) * vertex1;
00704     vertex2 = (transform) * vertex2;
00705     vertex3 = (transform) * vertex3;
00706     vertex4 = (transform) * vertex4;
00707     vertex5 = (transform) * vertex5;
00708     vertex6 = (transform) * vertex6;
00709     vertex7 = (transform) * vertex7;
00710     vertex8 = (transform) * vertex8;
00711 
00712     HepRepInstance* instance = getGeometryOrEventInstance(getCalHitType());
00713 
00714     addAttributes(instance, getCalHitType());
00715 
00716     setAttribute(instance, "DrawAs", G4String("Prism"));
00717         
00718     setVisibility(instance, trd);
00719     setLine(instance, trd);
00720     setColor(instance, getColorFor(trd));
00721 
00722     factory->createHepRepPoint(instance, vertex1.x(), vertex1.y(), vertex1.z());
00723     factory->createHepRepPoint(instance, vertex2.x(), vertex2.y(), vertex2.z());
00724     factory->createHepRepPoint(instance, vertex3.x(), vertex3.y(), vertex3.z());
00725     factory->createHepRepPoint(instance, vertex4.x(), vertex4.y(), vertex4.z());
00726     factory->createHepRepPoint(instance, vertex5.x(), vertex5.y(), vertex5.z());
00727     factory->createHepRepPoint(instance, vertex6.x(), vertex6.y(), vertex6.z());
00728     factory->createHepRepPoint(instance, vertex7.x(), vertex7.y(), vertex7.z());
00729     factory->createHepRepPoint(instance, vertex8.x(), vertex8.y(), vertex8.z());
00730 }
00731 
00732 void G4HepRepSceneHandler::AddSolid (const G4Trap& trap) { 
00733     if (dontWrite()) return;
00734     G4VSceneHandler::AddSolid (trap); 
00735 }
00736 
00737 void G4HepRepSceneHandler::AddSolid (const G4Sphere& sphere) { 
00738     if (dontWrite()) return;
00739     G4VSceneHandler::AddSolid (sphere); 
00740 }
00741 
00742 void G4HepRepSceneHandler::AddSolid (const G4Para& para) { 
00743     if (dontWrite()) return;
00744     G4VSceneHandler::AddSolid (para); 
00745 }
00746 
00747 void G4HepRepSceneHandler::AddSolid (const G4Torus& torus) { 
00748     if (dontWrite()) return;
00749     G4VSceneHandler::AddSolid (torus); 
00750 }
00751 
00752 void G4HepRepSceneHandler::AddSolid (const G4Polycone& polycone) { 
00753     if (dontWrite()) return;
00754     G4VSceneHandler::AddSolid (polycone); 
00755 }
00756 
00757 void G4HepRepSceneHandler::AddSolid (const G4Polyhedra& polyhedra) { 
00758     if (dontWrite()) return;
00759     G4VSceneHandler::AddSolid (polyhedra); 
00760 }
00761 
00762 void G4HepRepSceneHandler::AddSolid (const G4VSolid& solid) { 
00763     if (dontWrite()) return;
00764     G4VSceneHandler::AddSolid(solid); 
00765 }
00766 
00767 
00768 void G4HepRepSceneHandler::AddPrimitive (const G4Polyline& line) {
00769 
00770 #ifdef PDEBUG
00771     cout << "G4HepRepSceneHandler::AddPrimitive(G4Polyline&) " << line.size() << endl;
00772 #endif
00773     if (dontWrite()) return;
00774 
00775     if (fProcessing2D) {
00776       static G4bool warned = false;
00777       if (!warned) {
00778         warned = true;
00779         G4Exception
00780           ("G4HepRepSceneHandler::AddPrimitive (const G4Polyline&)",
00781            "vis-HepRep1001", JustWarning,
00782            "2D polylines not implemented.  Ignored.");
00783       }
00784       return;
00785     }
00786 
00787     HepRepInstance* instance = factory->createHepRepInstance(getEventInstance(), getTrajectoryType());
00788 
00789     addAttributes(instance, getTrajectoryType());
00790 
00791     setColor(instance, GetColor(line));
00792 
00793     setVisibility(instance, line);
00794 
00795     setLine(instance, line);
00796 
00797     for (size_t i=0; i < line.size(); i++) {
00798         G4Point3D vertex = transform * line[i];
00799         factory->createHepRepPoint(instance, vertex.x(), vertex.y(), vertex.z());
00800     }
00801 }
00802 
00803 
00804 void G4HepRepSceneHandler::AddPrimitive (const G4Polymarker& line) {
00805 
00806 #ifdef PDEBUG
00807     cout << "G4HepRepSceneHandler::AddPrimitive(G4Polymarker&) " << line.size() << endl;
00808 #endif
00809     if (dontWrite()) return;
00810 
00811     if (fProcessing2D) {
00812       static G4bool warned = false;
00813       if (!warned) {
00814         warned = true;
00815         G4Exception
00816           ("G4HepRepSceneHandler::AddPrimitive (const G4Polymarker&)",
00817            "vis-HepRep1002", JustWarning,
00818            "2D polymarkers not implemented.  Ignored.");
00819       }
00820       return;
00821     }
00822 
00823     HepRepInstance* instance = factory->createHepRepInstance(getEventInstance(), getHitType());
00824 
00825     addAttributes(instance, getHitType());
00826 
00827     setColor(instance, GetColor(line));
00828 
00829     setVisibility(instance, line);
00830 
00831     setMarker(instance, line);
00832 
00833     // Default MarkName is set to Circle for this Type.
00834     int mtype = line.GetMarkerType();
00835     
00836     // Cannot be case statement since line.xxx is not a constant
00837     if (mtype == line.dots) {
00838         setAttribute(instance, "Fill", true);
00839         setColor(instance, GetColor(line), G4String("FillColor"));
00840     } else if (mtype == line.circles) {
00841     } else if (line.squares) {
00842         setAttribute(instance, "MarkName", G4String("Box"));
00843     } else {
00844         // line.line + default
00845         setAttribute(instance, "MarkName", G4String("Plus"));
00846     }
00847 
00848     for (size_t i=0; i < line.size(); i++) {
00849         G4Point3D vertex = transform * line[i];
00850         factory->createHepRepPoint(instance, vertex.x(), vertex.y(), vertex.z());
00851     }
00852 }
00853 
00854 
00855 void G4HepRepSceneHandler::AddPrimitive (const G4Circle& circle) {
00856 #ifdef PDEBUG
00857     cout << "G4HepRepSceneHandler::AddPrimitive(G4Circle&) " << endl;
00858 #endif
00859     if (dontWrite()) return;
00860 
00861     if (fProcessing2D) {
00862       static G4bool warned = false;
00863       if (!warned) {
00864         warned = true;
00865         G4Exception
00866           ("G4HepRepSceneHandler::AddPrimitive (const G4Circle&)",
00867            "vis-HepRep1003", JustWarning,
00868            "2D circles not implemented.  Ignored.");
00869       }
00870       return;
00871     }
00872 
00873     HepRepInstance* instance = factory->createHepRepInstance(getEventInstance(), getHitType());
00874 
00875     addAttributes(instance, getHitType());
00876 
00877     G4Point3D center = transform * circle.GetPosition();
00878 
00879     setColor (instance, GetColor(circle));
00880 
00881     setVisibility(instance, circle);
00882 
00883     setMarker(instance, circle);
00884 
00885     factory->createHepRepPoint(instance, center.x(), center.y(), center.z());
00886 }
00887 
00888 
00889 void G4HepRepSceneHandler::AddPrimitive (const G4Polyhedron& polyhedron) {
00890 
00891 #ifdef PDEBUG
00892     cout << "G4HepRepSceneHandler::AddPrimitive(G4Polyhedron&) " << endl;
00893 #endif
00894     if (dontWrite()) return;
00895 
00896     if (fProcessing2D) {
00897       static G4bool warned = false;
00898       if (!warned) {
00899         warned = true;
00900         G4Exception
00901           ("G4HepRepSceneHandler::AddPrimitive (const G4Polyhedron&)",
00902            "vis-HepRep1004", JustWarning,
00903            "2D polyhedra not implemented.  Ignored.");
00904       }
00905       return;
00906     }
00907 
00908     G4Normal3D surfaceNormal;
00909     G4Point3D vertex;
00910 
00911     if (polyhedron.GetNoFacets()==0) return;
00912 
00913     HepRepInstance* instance = getGeometryOrEventInstance(getCalHitType());
00914         
00915     addAttributes(instance, getCalHitType());
00916 
00917     setVisibility(instance, polyhedron);
00918         
00919     G4int currentDepth = 0;
00920     G4PhysicalVolumeModel* pPVModel =
00921       dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
00922     if (pPVModel) currentDepth = pPVModel->GetCurrentDepth();
00923 
00924     G4bool notLastFace;
00925     do {
00926         HepRepInstance* face;
00927         if (isEventData()) {
00928             face = factory->createHepRepInstance(instance, getCalHitFaceType());
00929         } else {
00930             face = getGeometryInstance("*Face", currentDepth+1);
00931             setAttribute(face, "PickParent", true);
00932             setAttribute(face, "DrawAs", G4String("Polygon"));
00933         }
00934         
00935         setLine(face, polyhedron);
00936         setColor(face, GetColor(polyhedron));
00937         if (isEventData()) setColor(face, GetColor(polyhedron), G4String("FillColor"));
00938 
00939         notLastFace = polyhedron.GetNextNormal (surfaceNormal);
00940 
00941         G4int edgeFlag = 1;
00942         G4bool notLastEdge;
00943         do {
00944                 notLastEdge = polyhedron.GetNextVertex (vertex, edgeFlag);
00945             vertex = transform * vertex;
00946             factory->createHepRepPoint(face, vertex.x(), vertex.y(), vertex.z());
00947         } while (notLastEdge);
00948     } while (notLastFace);
00949 }
00950 
00951 
00952 void G4HepRepSceneHandler::AddPrimitive (const G4Text&) {
00953 #ifdef PDEBUG
00954     cout << "G4HepRepSceneHandler::AddPrimitive(G4Text&) " << endl;
00955 #endif
00956     if (dontWrite()) return;
00957 
00958     /*** You may need this
00959     if (fProcessing2D) {
00960       static G4bool warned = false;
00961       if (!warned) {
00962         warned = true;
00963         G4Exception
00964           ("G4HepRepSceneHandler::AddPrimitive (const G4Text&)",
00965            "vis-HepRep1005", JustWarning,
00966            "2D text not implemented.  Ignored.");
00967       }
00968       return;
00969     }
00970     ***/
00971 
00972     cout << "G4HepRepSceneHandler::AddPrimitive G4Text : not yet implemented. " << endl;
00973 }
00974 
00975 
00976 void G4HepRepSceneHandler::AddPrimitive (const G4Square& square) {
00977 #ifdef PDEBUG
00978     cout << "G4HepRepSceneHandler::AddPrimitive(G4Square&) " << endl;
00979 #endif
00980     if (dontWrite()) return;
00981 
00982     if (fProcessing2D) {
00983       static G4bool warned = false;
00984       if (!warned) {
00985         warned = true;
00986         G4Exception
00987           ("G4HepRepSceneHandler::AddPrimitive (const G4Square&)",
00988            "vis-HepRep1006", JustWarning,
00989            "2D squares not implemented.  Ignored.");
00990       }
00991       return;
00992     }
00993 
00994     HepRepInstance* instance = factory->createHepRepInstance(getEventInstance(), getHitType());
00995 
00996     addAttributes(instance, getHitType());
00997 
00998     G4Point3D center = transform * square.GetPosition();
00999 
01000     setColor (instance, getColorFor(square));
01001 
01002     setVisibility(instance, square);
01003 
01004     setMarker(instance, square);
01005 
01006     factory->createHepRepPoint(instance, center.x(), center.y(), center.z());
01007 }
01008 
01009 void G4HepRepSceneHandler::AddPrimitive (const G4NURBS&) {
01010 #ifdef PDEBUG
01011     cout << "G4HepRepSceneHandler::AddPrimitive(G4NURBS&) " << endl;
01012 #endif
01013     if (dontWrite()) return;
01014 
01015     /*** You may need this
01016     if (fProcessing2D) {
01017       static G4bool warned = false;
01018       if (!warned) {
01019         warned = true;
01020         G4Exception
01021           ("G4HepRepSceneHandler::AddPrimitive (const G4NURBS&)",
01022            "vis-HepRep1007", JustWarning,
01023            "2D NURBS not implemented.  Ignored.");
01024       }
01025       return;
01026     }
01027     ***/
01028 
01029     cout << "G4HepRepSceneHandler::AddPrimitive G4NURBS : not yet implemented. " << endl;
01030 }
01031 
01032 void G4HepRepSceneHandler::AddPrimitive (const G4Scale& scale) { 
01033     if (dontWrite()) return;
01034     G4VSceneHandler::AddPrimitive(scale); 
01035 }
01036 
01037 void G4HepRepSceneHandler::AddCompound (const G4VTrajectory& trajectory) {
01038 #ifdef PDEBUG
01039     cout << "G4HepRepSceneHandler::AddCompound(G4VTrajectory&) " << endl;
01040 #endif
01041     if (dontWrite()) return;
01042     
01043     currentTrack = &trajectory;
01044     G4VSceneHandler::AddCompound(trajectory); 
01045     currentTrack = NULL;    
01046 }
01047 
01048 
01049 void G4HepRepSceneHandler::AddCompound (const G4VHit& hit) { 
01050 #ifdef PDEBUG
01051     cout << "G4HepRepSceneHandler::AddCompound(G4VHit&) " << endl;
01052 #endif
01053     if (dontWrite()) return;
01054     
01055     currentHit = &hit;
01056     G4VSceneHandler::AddCompound(hit); 
01057     currentHit = NULL;
01058 }
01059 
01060 void G4HepRepSceneHandler::PreAddSolid (const G4Transform3D& objectTransformation,
01061                                   const G4VisAttributes& visAttribs) {
01062 
01063     G4VSceneHandler::PreAddSolid (objectTransformation, visAttribs);
01064 
01065     transform = objectTransformation;
01066 #ifdef SDEBUG
01067     cout << "G4HepRepSceneHandler::PreAddSolid(G4Transform3D&, G4VisAttributes&)" << endl;
01068 #endif
01069 }
01070 
01071 
01072 void G4HepRepSceneHandler::PostAddSolid () {
01073 #ifdef SDEBUG
01074     cout << "G4HepRepSceneHandler::PostAddSolid()" << endl;
01075 #endif
01076     G4VSceneHandler::PostAddSolid();
01077 }
01078 
01079 
01080 void G4HepRepSceneHandler::BeginPrimitives (const G4Transform3D& objectTransformation) {
01081 #ifdef SDEBUG
01082     cout << "G4HepRepSceneHandler::BeginPrimitives(G4Transform3D&)" << endl;
01083 #endif
01084 
01085     G4VSceneHandler::BeginPrimitives (objectTransformation);
01086     transform = objectTransformation;
01087 }
01088 
01089 
01090 void G4HepRepSceneHandler::EndPrimitives () {
01091 #ifdef SDEBUG
01092     cout << "G4HepRepSceneHandler::EndPrimitives" << endl;
01093 #endif
01094     G4VSceneHandler::EndPrimitives ();
01095 }
01096 
01097 
01098 G4bool G4HepRepSceneHandler::dontWrite() {      
01099         G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01100     return !( messenger->writeInvisibles() || (fpVisAttribs ? (bool)fpVisAttribs->IsVisible() : true));
01101 }
01102 
01103 void G4HepRepSceneHandler::setColor (HepRepAttribute *attribute,
01104                                       const G4Color& color,
01105                                       const G4String& key) {
01106 #ifdef CDEBUG
01107     cout << "G4HepRepSceneHandler::setColor : red : " << color.GetRed ()   <<
01108                                   " green : " << color.GetGreen () <<
01109                                   " blue : " << color.GetBlue ()   << endl;
01110 #endif
01111 
01112     setAttribute(attribute, key, color.GetRed(), color.GetGreen(), color.GetBlue(), color.GetAlpha());
01113 }
01114 
01115 G4Color G4HepRepSceneHandler::getColorFor (const G4VSolid& /* solid */) {
01116     return fpVisAttribs ? fpVisAttribs->GetColor() : GetColor(NULL);    
01117 }
01118 
01119 G4Color G4HepRepSceneHandler::getColorFor (const G4Visible& visible) {
01120     return GetColor(visible);
01121 }
01122 
01123 void G4HepRepSceneHandler::setVisibility (HepRepAttribute *attribute, const G4VSolid& /* solid */) {
01124     setAttribute(attribute, "Visibility", (fpVisAttribs ? (bool)fpVisAttribs->IsVisible() : true));
01125 }
01126 
01127 void G4HepRepSceneHandler::setVisibility ( HepRepAttribute *attribute, const G4Visible& visible) {
01128     const G4VisAttributes* atts = visible.GetVisAttributes();
01129 
01130     setAttribute(attribute, "Visibility", (atts && (atts->IsVisible()==0)) ? false : true);
01131 }
01132 
01133 void G4HepRepSceneHandler::setLine (HepRepAttribute *attribute, const G4VSolid& /* solid*/) {
01134     setAttribute(attribute, "LineWidth", 1.0);
01135 }
01136 
01137 void G4HepRepSceneHandler::setLine (HepRepAttribute *attribute, const G4Visible& visible) {
01138     const G4VisAttributes* atts = visible.GetVisAttributes();
01139     
01140     setAttribute(attribute, "LineWidth", (atts != NULL) ? atts->GetLineWidth() : 1.0);
01141     
01142     if (atts != NULL) {
01143         switch (atts->GetLineStyle()) {
01144             case G4VisAttributes::dotted:
01145                 setAttribute(attribute, "LineStyle", G4String("Dotted"));
01146                 break;
01147             case G4VisAttributes::dashed:
01148                 setAttribute(attribute, "LineStyle", G4String("Dashed"));
01149                 break;
01150             case G4VisAttributes::unbroken:
01151             default:
01152                 break;
01153         }
01154     }
01155 }
01156 
01157 void G4HepRepSceneHandler::setMarker (HepRepAttribute *attribute, const G4VMarker& marker) {
01158     MarkerSizeType markerType;
01159     G4double size = GetMarkerRadius( marker , markerType );
01160 
01161     setAttribute(attribute, "MarkSize", size);
01162 
01163     if (markerType == screen) setAttribute(attribute, "MarkType", G4String("Symbol"));
01164     if (marker.GetFillStyle() == G4VMarker::noFill) {
01165         setAttribute(attribute, "Fill", false);
01166     } else {
01167         setColor(attribute, GetColor(marker), G4String("FillColor"));
01168     }
01169 }
01170 
01171 void G4HepRepSceneHandler::addAttributes(HepRepInstance* instance, HepRepType* type) {
01172     if (currentHit) {
01173         vector<G4AttValue>* hitAttValues = currentHit->CreateAttValues();
01174         const map<G4String,G4AttDef>* hitAttDefs = currentHit->GetAttDefs();
01175 
01176         addAttDefs(getHitType(), hitAttDefs);
01177 
01178         // these attValues are non-standard, so can only be added when we have the attDef.
01179         type->addAttValue("LVol", G4String(""));
01180         type->addAttValue("HitType", G4String(""));
01181         type->addAttValue("ID", -1);
01182         type->addAttValue("Column", -1);
01183         type->addAttValue("Row", -1);
01184         type->addAttValue("Energy", 0.0);
01185         type->addAttValue("Pos", G4String(""));
01186 
01187         addAttVals(instance, hitAttDefs, hitAttValues);
01188     
01189         delete hitAttValues;
01190     
01191     } else if (currentTrack) {
01192         vector<G4AttValue>* trajectoryAttValues = currentTrack->CreateAttValues();
01193         const map<G4String,G4AttDef>* trajectoryAttDefs = currentTrack->GetAttDefs();
01194 
01195         addAttDefs(type, trajectoryAttDefs);
01196     
01197         // these attValues are non-standard, so can only be added when we have the attDef.
01198         type->addAttValue("Ch", 0.0);
01199         type->addAttValue("Color", 1.0, 1.0, 1.0, 1.0);
01200         type->addAttValue("ID", -1);
01201         type->addAttValue("IMom", G4String(""));
01202         type->addAttValue("IMag", 0.0);
01203         type->addAttValue("PDG", -1);
01204         type->addAttValue("PN", G4String(""));
01205         type->addAttValue("PID", -1);
01206 
01207         addAttVals(instance, trajectoryAttDefs, trajectoryAttValues);
01208     
01209         delete trajectoryAttValues;
01210         
01211     }
01212 }
01213 
01214 void G4HepRepSceneHandler::setAttribute(HepRepAttribute* attribute, G4String name, G4String value) {
01215     HepRepAttValue* attValue = attribute->getAttValue(name);
01216     if ((attValue == NULL) || (attValue->getString() != value)) {
01217         HepRepPoint* point = dynamic_cast<HepRepPoint*>(attribute);
01218         if (point != NULL) {
01219             if (point->getInstance()->getAttValueFromNode(name) == NULL) {
01220                 attribute = point->getInstance();
01221             }
01222         }
01223         
01224         HepRepInstance* instance = dynamic_cast<HepRepInstance*>(attribute);
01225         if (instance != NULL) {
01226             // look for definition on type (node only)
01227             if (instance->getType()->getAttValueFromNode(name) == NULL) {
01228                 attribute = instance->getType();
01229             }   
01230         }
01231         
01232         attribute->addAttValue(name, value);
01233     }
01234 }
01235 
01236 void G4HepRepSceneHandler::setAttribute(HepRepAttribute* attribute, G4String name, bool value) {
01237     HepRepAttValue* attValue = attribute->getAttValue(name);
01238     if ((attValue == NULL) || (attValue->getBoolean() != value)) {
01239         HepRepPoint* point = dynamic_cast<HepRepPoint*>(attribute);
01240         if (point != NULL) {
01241             if (point->getInstance()->getAttValueFromNode(name) == NULL) {
01242                 attribute = point->getInstance();
01243             }
01244         }
01245         
01246         HepRepInstance* instance = dynamic_cast<HepRepInstance*>(attribute);
01247         if (instance != NULL) {
01248             // look for definition on type (node only)
01249             if (instance->getType()->getAttValueFromNode(name) == NULL) {
01250                 attribute = instance->getType();
01251             }   
01252         }
01253         
01254         attribute->addAttValue(name, value);
01255     }
01256 }
01257 
01258 void G4HepRepSceneHandler::setAttribute(HepRepAttribute* attribute, G4String name, double value) {
01259     HepRepAttValue* attValue = attribute->getAttValue(name);
01260     if ((attValue == NULL) || (attValue->getDouble() != value)) {
01261         HepRepPoint* point = dynamic_cast<HepRepPoint*>(attribute);
01262         if (point != NULL) {
01263             if (point->getInstance()->getAttValueFromNode(name) == NULL) {
01264                 attribute = point->getInstance();
01265             }
01266         }
01267         
01268         HepRepInstance* instance = dynamic_cast<HepRepInstance*>(attribute);
01269         if (instance != NULL) {
01270             // look for definition on type (node only)
01271             if (instance->getType()->getAttValueFromNode(name) == NULL) {
01272                 attribute = instance->getType();
01273             }   
01274         }
01275         
01276         attribute->addAttValue(name, value);
01277     }
01278 }
01279 
01280 void G4HepRepSceneHandler::setAttribute(HepRepAttribute* attribute, G4String name, int value) {
01281     HepRepAttValue* attValue = attribute->getAttValue(name);
01282     if ((attValue == NULL) || (attValue->getInteger() != value)) {
01283         HepRepPoint* point = dynamic_cast<HepRepPoint*>(attribute);
01284         if (point != NULL) {
01285             if (point->getInstance()->getAttValueFromNode(name) == NULL) {
01286                 attribute = point->getInstance();
01287             }
01288         }
01289         
01290         HepRepInstance* instance = dynamic_cast<HepRepInstance*>(attribute);
01291         if (instance != NULL) {
01292             // look for definition on type (node only)
01293             if (instance->getType()->getAttValueFromNode(name) == NULL) {
01294                 attribute = instance->getType();
01295             }   
01296         }
01297         
01298         attribute->addAttValue(name, value);
01299     }
01300 }
01301 
01302 void G4HepRepSceneHandler::setAttribute(HepRepAttribute* attribute, G4String name, double red, double green, double blue, double alpha) {
01303     HepRepAttValue* attValue = attribute->getAttValue(name);
01304     vector<double> color;
01305     if (attValue != NULL) color = attValue->getColor();
01306     if ((color.size() == 0) ||
01307         (color[0] != red) ||
01308         (color[1] != green) ||
01309         (color[2] != blue) ||
01310         ((color.size() > 3) && (color[3] != alpha))) {
01311         
01312         HepRepPoint* point = dynamic_cast<HepRepPoint*>(attribute);
01313         if (point != NULL) {
01314             if (point->getInstance()->getAttValueFromNode(name) == NULL) {
01315                 attribute = point->getInstance();
01316             }
01317         }
01318         
01319         HepRepInstance* instance = dynamic_cast<HepRepInstance*>(attribute);
01320         if (instance != NULL) {
01321             // look for definition on type (node only)
01322             if (instance->getType()->getAttValueFromNode(name) == NULL) {
01323                 attribute = instance->getType();
01324             }   
01325         }
01326         
01327         attribute->addAttValue(name, red, green, blue, alpha);
01328     }
01329 }
01330 
01331 void G4HepRepSceneHandler::addAttDefs(HepRepDefinition* definition, const map<G4String,G4AttDef>* attDefs) {
01332     if (attDefs == NULL) return;
01333 
01334     // Specify additional attribute definitions.
01335     map<G4String,G4AttDef>::const_iterator attDefIterator = attDefs->begin();
01336     while (attDefIterator != attDefs->end()) {
01337         definition->addAttDef(attDefIterator->first, attDefIterator->second.GetDesc(),
01338                         attDefIterator->second.GetCategory(), attDefIterator->second.GetExtra());
01339         attDefIterator++;
01340     }
01341 }
01342 
01343 void G4HepRepSceneHandler::addAttVals(HepRepAttribute* attribute, const map<G4String,G4AttDef>* attDefs, vector<G4AttValue>* attValues) {
01344     if (attValues == NULL) return;
01345 
01346     // Copy the instance's G4AttValues to HepRepAttValues.
01347     for (vector<G4AttValue>::iterator attValIterator = attValues->begin(); attValIterator != attValues->end(); attValIterator++) {        
01348         G4String name = attValIterator->GetName();
01349 
01350         HepRepPoint* point = dynamic_cast<HepRepPoint*>(attribute);
01351         if ((name == "Pos") && (point != NULL)) {
01352             G4String pos = attValIterator->GetValue();
01353 //            cout << "Pos* " << pos << endl;
01354             int is = 0;
01355             int in = 0;
01356             int im = 0;
01357             G4String unit;
01358             for (unsigned int i=0; i<pos.length(); i++) {
01359                 if (pos[i] == ' ') {
01360                     if (in == 0) {
01361                         // first coordinate
01362                         double factor = atof(pos.substr(is, i-is).c_str())/point->getX();
01363                         im = (int)(std::log10(factor)+((factor < 1) ? -0.5 : 0.5));
01364 //                        cout << factor << ", " << im << endl;
01365                     } else if (in == 3) {
01366                         // unit
01367                         unit = pos.substr(is, i-is);
01368                         if (unit == G4String("mum")) {
01369                             im += -6;
01370                         } else if (unit == G4String("mm")) {
01371                             im += -3;
01372                         } else if (unit == G4String("cm")) {
01373                             im += -2;
01374                         } else if (unit == G4String("m")) {
01375                             im += 0;
01376                         } else if (unit == G4String("km")) {
01377                             im += 3;
01378                         } else {
01379                             cerr << "HepRepSceneHandler: Unrecognized Unit: '" << unit << "'" << endl;
01380                         }
01381                     }
01382                     is = i+1;
01383                     in++;
01384                 }
01385             }
01386             switch(im) {
01387                 case -6:
01388                     unit = G4String("mum");
01389                     break;
01390                 case -3:
01391                     unit = G4String("mm");
01392                     break;
01393                 case -2:
01394                     unit = G4String("cm");
01395                     break;
01396                 case 0:
01397                     unit = G4String("m");
01398                     break;
01399                 case 3:
01400                     unit = G4String("km");
01401                     break;
01402                 default:
01403                     cerr << "HepRepSceneHandler: No valid unit found for im: " << im << endl;
01404                     unit = G4String("*im");                        
01405                     break;
01406             }
01407 //            cout << "U: " << unit << endl;
01408             setAttribute(attribute, G4String("PointUnit"), unit);
01409             continue;
01410         }
01411                 
01412         // NTP already in points being written
01413         if (name == "NTP") continue;
01414         
01415         // find type of attribute using def
01416         const map<G4String,G4AttDef>::const_iterator attDefIterator = attDefs->find(name);
01417         G4String type = attDefIterator->second.GetValueType();
01418                 
01419         // set based on type
01420         if ((type == "G4double") || (type == "double")) {   
01421             setAttribute(attribute, attValIterator->GetName(), atof(attValIterator->GetValue()));           
01422         } else if ((type == "G4int") || (type == "int")) {  
01423             setAttribute(attribute, attValIterator->GetName(), atoi(attValIterator->GetValue()));           
01424         } else { // G4String, string and others
01425             setAttribute(attribute, attValIterator->GetName(), attValIterator->GetValue());
01426         }
01427     }
01428 }
01429 
01430 
01431 bool G4HepRepSceneHandler::isEventData () {
01432     G4PhysicalVolumeModel* pPVModel =
01433       dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
01434     return !pPVModel || fReadyForTransients || currentHit || currentTrack;
01435 }
01436 
01437 void G4HepRepSceneHandler::addTopLevelAttributes(HepRepType* type) {
01438     
01439     // Some non-standard attributes
01440     type->addAttDef(  "Generator", "Generator of the file", "General", "");
01441     type->addAttValue("Generator", G4String("Geant4"));
01442 
01443     type->addAttDef(  "GeneratorVersion", "Version of the Generator", "General", "");
01444     G4String versionString = G4Version;
01445     versionString = versionString.substr(1,versionString.size()-2);
01446     versionString = " Geant4 version " + versionString + "   " + G4Date;
01447     type->addAttValue("GeneratorVersion", versionString);
01448     
01449     const G4ViewParameters parameters = GetCurrentViewer()->GetViewParameters();
01450     const G4Vector3D& viewPointDirection = parameters.GetViewpointDirection();    
01451     type->addAttDef(  "ViewTheta", "Theta of initial suggested viewpoint", "Draw", "rad");
01452     type->addAttValue("ViewTheta", viewPointDirection.theta());    
01453     
01454     type->addAttDef(  "ViewPhi", "Phi of initial suggested viewpoint", "Draw", "rad");
01455     type->addAttValue("ViewPhi", viewPointDirection.phi());    
01456     
01457     type->addAttDef(  "ViewScale", "Scale of initial suggested viewpoint", "Draw", "");
01458     type->addAttValue("ViewScale", parameters.GetZoomFactor());    
01459     
01460 // FIXME, no way to set these
01461     type->addAttDef(  "ViewTranslateX", "Translate in X of initial suggested viewpoint", "Draw", "");
01462     type->addAttValue("ViewTranslateX", 0.0);    
01463     
01464     type->addAttDef(  "ViewTranslateY", "Translate in Y of initial suggested viewpoint", "Draw", "");
01465     type->addAttValue("ViewTranslateY", 0.0);    
01466     
01467     type->addAttDef(  "ViewTranslateZ", "Translate in Z of initial suggested viewpoint", "Draw", "");
01468     type->addAttValue("ViewTranslateZ", 0.0);    
01469     
01470     type->addAttDef(  "PointUnit", "Length", "Physics", "");
01471     type->addAttValue("PointUnit", G4String("m"));
01472         
01473         G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01474 
01475     type->addAttDef(  "UseSolids", "Use HepRep Solids rather than Geant4 Primitives", "Draw", "");
01476     type->addAttValue("UseSolids", messenger->useSolids());
01477 
01478     type->addAttDef(  "WriteInvisibles", "Write Invisible Objects", "Draw", "");
01479     type->addAttValue("WriteInvisibles", messenger->writeInvisibles());
01480 }           
01481 
01482 
01483 HepRepInstance* G4HepRepSceneHandler::getGeometryOrEventInstance(HepRepType* type) {
01484     if (isEventData()) {
01485       return factory->createHepRepInstance(getEventInstance(), type);
01486     } else {
01487       G4PhysicalVolumeModel* pPVModel =
01488         dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
01489       G4LogicalVolume* pCurrentLV = pPVModel->GetCurrentLV();
01490       G4int currentDepth = pPVModel->GetCurrentDepth();
01491       G4Material* pCurrentMaterial = pPVModel->GetCurrentMaterial();
01492       return getGeometryInstance(pCurrentLV, pCurrentMaterial, currentDepth);
01493     }
01494 }
01495 
01496 HepRep* G4HepRepSceneHandler::getHepRep() {
01497     if (_heprep == NULL) {
01498         // Create the HepRep that holds the Trees.
01499         _heprep = factory->createHepRep();
01500     }
01501     return _heprep;
01502 }   
01503 
01504 HepRep* G4HepRepSceneHandler::getHepRepGeometry() {
01505     if (_heprepGeometry == NULL) {
01506         // Create the HepRep that holds the Trees.
01507         _heprepGeometry = factory->createHepRep();
01508     }
01509     return _heprepGeometry;
01510 }   
01511 
01512 HepRepInstanceTree* G4HepRepSceneHandler::getGeometryInstanceTree() {
01513     if (_geometryInstanceTree == NULL) {
01514         // Create the Geometry InstanceTree.
01515         _geometryInstanceTree = factory->createHepRepInstanceTree("G4GeometryData", "1.0", getGeometryTypeTree());
01516                 
01517                 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01518         if ( messenger->appendGeometry()) {
01519             getHepRep()->addInstanceTree(_geometryInstanceTree);
01520         } else {
01521             getHepRepGeometry()->addInstanceTree(_geometryInstanceTree);
01522         }
01523     }
01524     return _geometryInstanceTree;
01525 }
01526 
01527 HepRepInstance* G4HepRepSceneHandler::getGeometryRootInstance() {
01528     if (_geometryRootInstance == NULL) {
01529         // Create the top level Geometry Instance.
01530         _geometryRootInstance = factory->createHepRepInstance(getGeometryInstanceTree(), getGeometryRootType());
01531     }
01532     return _geometryRootInstance;
01533 }
01534 
01535 HepRepInstance* G4HepRepSceneHandler::getGeometryInstance(G4LogicalVolume* volume, G4Material* material, int depth) {
01536     HepRepInstance* instance = getGeometryInstance(volume->GetName(), depth);
01537 
01538     setAttribute(instance, "LVol",       volume->GetName());
01539     G4Region* region = volume->GetRegion();
01540     G4String regionName = region? region->GetName(): G4String("No region");
01541     setAttribute(instance, "Region",     regionName);
01542     setAttribute(instance, "RootRegion", volume->IsRootRegion());
01543     setAttribute(instance, "Solid",      volume->GetSolid()->GetName());
01544     setAttribute(instance, "EType",      volume->GetSolid()->GetEntityType());
01545     G4String matName = material? material->GetName(): G4String("No material");
01546     setAttribute(instance, "Material",   matName );
01547     G4double matDensity = material? material->GetDensity(): 0.;
01548     setAttribute(instance, "Density",    matDensity);
01549     G4double matRadlen = material? material->GetRadlen(): 0.;
01550     setAttribute(instance, "Radlen",     matRadlen);
01551     
01552     G4State matState = material? material->GetState(): kStateUndefined;
01553     G4String state = materialState[matState];
01554     setAttribute(instance, "State", state);
01555     
01556     return instance;
01557 }
01558 
01559 HepRepInstance* G4HepRepSceneHandler::getGeometryInstance(G4String volumeName, int depth) {
01560     // no extra checks since these are done in the geometryType already
01561 
01562     // adjust depth, also pop the current instance
01563     while ((int)_geometryInstance.size() > depth) {
01564         _geometryInstance.pop_back();
01565     }
01566     
01567     // get parent
01568     HepRepInstance* parent = (_geometryInstance.empty()) ? getGeometryRootInstance() : _geometryInstance.back();
01569     
01570     // get type
01571     HepRepType* type = getGeometryType(volumeName, depth);
01572     
01573     // create instance
01574     HepRepInstance* instance = factory->createHepRepInstance(parent, type);
01575     _geometryInstance.push_back(instance);
01576 
01577     return instance;
01578 }
01579 
01580 HepRepTypeTree* G4HepRepSceneHandler::getGeometryTypeTree() {
01581     if (_geometryTypeTree == NULL) {
01582         // Create the Geometry TypeTree.
01583         HepRepTreeID* geometryTreeID = factory->createHepRepTreeID("G4GeometryTypes", "1.0");
01584         _geometryTypeTree = factory->createHepRepTypeTree(geometryTreeID);
01585                 
01586                 G4HepRepMessenger* messenger = G4HepRepMessenger::GetInstance();
01587         if ( messenger->appendGeometry()) {
01588             getHepRep()->addTypeTree(_geometryTypeTree);
01589         } else {
01590             getHepRepGeometry()->addTypeTree(_geometryTypeTree);
01591         }
01592     }
01593     return _geometryTypeTree;
01594 }
01595 
01596 HepRepType* G4HepRepSceneHandler::getGeometryRootType() {
01597     if (_geometryRootType == NULL) {
01598         // Create the top level Geometry Type.
01599         _geometryRootType = factory->createHepRepType(getGeometryTypeTree(), rootVolumeName);
01600         _geometryRootType->addAttValue("Layer", geometryLayer);
01601     
01602         // Add attdefs used by all geometry types.
01603         _geometryRootType->addAttDef  ("LVol", "Logical Volume", "Physics","");
01604         _geometryRootType->addAttValue("LVol", G4String(""));
01605         _geometryRootType->addAttDef  ("Region", "Cuts Region", "Physics","");
01606         _geometryRootType->addAttValue("Region", G4String(""));
01607         _geometryRootType->addAttDef  ("RootRegion", "Root Region", "Physics","");
01608         _geometryRootType->addAttValue("RootRegion", false);
01609         _geometryRootType->addAttDef  ("Solid", "Solid Name", "Physics","");
01610         _geometryRootType->addAttValue("Solid", G4String(""));
01611         _geometryRootType->addAttDef  ("EType", "Entity Type", "Physics","");
01612         _geometryRootType->addAttValue("EType", G4String("G4Box"));
01613         _geometryRootType->addAttDef  ("Material", "Material Name", "Physics","");
01614         _geometryRootType->addAttValue("Material", G4String("Air"));
01615         _geometryRootType->addAttDef  ("Density", "Material Density", "Physics","");
01616         _geometryRootType->addAttValue("Density", 0.0);
01617         _geometryRootType->addAttDef  ("State", "Material State", "Physics","");
01618         _geometryRootType->addAttValue("State", G4String("Gas"));
01619         _geometryRootType->addAttDef  ("Radlen", "Material Radiation Length", "Physics","");
01620         _geometryRootType->addAttValue("Radlen", 0.0);
01621         
01622         // add defaults for Geometry
01623         _geometryRootType->addAttValue("Color", 0.8, 0.8, 0.8, 1.0);
01624         _geometryRootType->addAttValue("Visibility", true);
01625         _geometryRootType->addAttValue("FillColor", 0.8, 0.8, 0.8, 1.0);
01626         _geometryRootType->addAttValue("LineWidth", 1.0);
01627         _geometryRootType->addAttValue("DrawAs", G4String("Polygon"));
01628         _geometryRootType->addAttValue("PickParent", false);
01629         _geometryRootType->addAttValue("ShowParentAttributes", true);
01630         
01631         _geometryRootType->addAttValue("MarkSizeMultiplier", 4.0);
01632         _geometryRootType->addAttValue("LineWidthMultiplier", 1.0);
01633 
01634         addTopLevelAttributes(_geometryRootType);       
01635                 
01636         _geometryType["/"+_geometryRootType->getName()] = _geometryRootType;
01637     }
01638     return _geometryRootType;
01639 }
01640 
01641 HepRepType* G4HepRepSceneHandler::getGeometryType(G4String volumeName, int depth) {
01642     // make sure we have a root
01643     getGeometryRootType();   
01644 
01645     // construct the full name for this volume
01646     G4String name = getFullTypeName(volumeName, depth);
01647     
01648     // lookup type and create if necessary
01649     HepRepType* type = _geometryType[name];    
01650     if (type == NULL) {
01651         G4String parentName = getParentTypeName(depth);
01652         HepRepType* parentType = _geometryType[parentName];
01653         // HepRep uses hierarchical names
01654         type = factory->createHepRepType(parentType, volumeName);
01655         _geometryType[name] = type;
01656     }
01657     return type;   
01658 }
01659 
01660 G4String G4HepRepSceneHandler::getFullTypeName(G4String volumeName, int depth) {
01661     // check for name depth
01662     if (depth > (int)_geometryTypeName.size()) {
01663         // there is a problem, book this type under problems
01664         G4String problem = "HierarchyProblem";
01665         if (_geometryType["/"+problem] == NULL) {
01666             // HepRep uses hierarchical names
01667             HepRepType* type = factory->createHepRepType(getGeometryRootType(), problem);
01668             _geometryType["/"+problem] = type;
01669         }
01670         return "/" + problem + "/" + volumeName;
01671     }
01672     
01673     // adjust name depth, also pop the current volumeName
01674     while ((int)_geometryTypeName.size() > depth) {
01675         _geometryTypeName.pop_back();
01676     }
01677     
01678     // construct full name and push it
01679     G4String name = (_geometryTypeName.empty()) ? G4String("/"+rootVolumeName) : _geometryTypeName.back();
01680     name = name + "/" + volumeName;
01681     _geometryTypeName.push_back(name);
01682     return name;
01683 }
01684 
01685 G4String G4HepRepSceneHandler::getParentTypeName(int depth) {
01686     return (depth >= 1) ? _geometryTypeName[depth-1] : G4String("/"+rootVolumeName);
01687 }       
01688 
01689 HepRepInstanceTree* G4HepRepSceneHandler::getEventInstanceTree() {
01690     if (_eventInstanceTree == NULL) {
01691         // Create the Event InstanceTree.
01692         _eventInstanceTree = factory->createHepRepInstanceTree("G4EventData", "1.0", getEventTypeTree());
01693         getHepRep()->addInstanceTree(_eventInstanceTree);
01694     }
01695     return _eventInstanceTree;
01696 }
01697 
01698 HepRepInstance* G4HepRepSceneHandler::getEventInstance() {
01699     if (_eventInstance == NULL) {
01700         // Create the top level Event Instance.
01701         _eventInstance = factory->createHepRepInstance(getEventInstanceTree(), getEventType());
01702     }
01703     return _eventInstance;
01704 }
01705 
01706 HepRepTypeTree* G4HepRepSceneHandler::getEventTypeTree() {
01707     if (_eventTypeTree == NULL) {
01708         // Create the Event TypeTree.
01709         HepRepTreeID* eventTreeID = factory->createHepRepTreeID("G4EventTypes", "1.0");
01710         _eventTypeTree = factory->createHepRepTypeTree(eventTreeID);
01711         getHepRep()->addTypeTree(_eventTypeTree);    
01712     }
01713     
01714     return _eventTypeTree;
01715 }
01716 
01717 HepRepType* G4HepRepSceneHandler::getEventType() {
01718     if (_eventType == NULL) {
01719         // Create the top level Event Type.
01720         _eventType = factory->createHepRepType(getEventTypeTree(), "Event");
01721         _eventType->addAttValue("Layer", eventLayer);
01722 
01723         // add defaults for Events
01724         _eventType->addAttValue("Visibility", true);
01725         _eventType->addAttValue("Color", 1.0, 1.0, 1.0, 1.0);
01726         _eventType->addAttValue("FillColor", 1.0, 1.0, 1.0, 1.0);
01727         _eventType->addAttValue("LineWidth", 1.0);
01728         _eventType->addAttValue("HasFrame", true);
01729         _eventType->addAttValue("PickParent", false);
01730         _eventType->addAttValue("ShowParentAttributes", false);
01731 
01732         _eventType->addAttValue("MarkSizeMultiplier", 4.0);
01733         _eventType->addAttValue("LineWidthMultiplier", 1.0);
01734 
01735         addTopLevelAttributes(_eventType);
01736     }
01737     
01738     return _eventType;
01739 }
01740 
01741 HepRepType* G4HepRepSceneHandler::getTrajectoryType() {
01742     if (_trajectoryType == NULL) {
01743         _trajectoryType = factory->createHepRepType(getEventType(), "Trajectory");
01744         
01745         _trajectoryType->addAttValue("Layer", trajectoryLayer);
01746         _trajectoryType->addAttValue("DrawAs", G4String("Line"));
01747 
01748         _trajectoryType->addAttValue("LineWidthMultiplier", 2.0);
01749         
01750         // attributes to draw the points of a track as markers.
01751         _trajectoryType->addAttValue("MarkName", G4String("Box"));
01752         _trajectoryType->addAttValue("MarkSize", 4);
01753         _trajectoryType->addAttValue("MarkType", G4String("Symbol"));
01754         _trajectoryType->addAttValue("Fill", true);
01755     }
01756     return _trajectoryType;
01757 }
01758 
01759 HepRepType* G4HepRepSceneHandler::getHitType() {
01760     if (_hitType == NULL) {
01761         _hitType = factory->createHepRepType(getEventType(), "Hit");
01762         _hitType->addAttValue("Layer", hitLayer);
01763         _hitType->addAttValue("DrawAs", G4String("Point"));
01764         _hitType->addAttValue("MarkName", G4String("Box"));
01765         _hitType->addAttValue("MarkSize", 4.0);
01766         _hitType->addAttValue("MarkType", G4String("Symbol"));
01767         _hitType->addAttValue("Fill", true);
01768     }
01769     return _hitType;
01770 }
01771 
01772 HepRepType* G4HepRepSceneHandler::getCalHitType() {
01773     if (_calHitType == NULL) {
01774         _calHitType = factory->createHepRepType(getEventType(), "CalHit");
01775         _calHitType->addAttValue("Layer", calHitLayer);
01776         _calHitType->addAttValue("Fill", true);
01777         _calHitType->addAttValue("DrawAs", G4String("Polygon"));
01778     }
01779     return _calHitType;
01780 }
01781 
01782 HepRepType* G4HepRepSceneHandler::getCalHitFaceType() {
01783     if (_calHitFaceType == NULL) {
01784         _calHitFaceType = factory->createHepRepType(getCalHitType(), "CalHitFace");
01785         _calHitFaceType->addAttValue("PickParent", true);
01786     }
01787     return _calHitFaceType;
01788 }
01789 

Generated on Mon May 27 17:48:30 2013 for Geant4 by  doxygen 1.4.7