source: trunk/source/visualization/gMocren/src/G4GMocrenFileSceneHandler.cc @ 1274

Last change on this file since 1274 was 1274, checked in by garnier, 14 years ago

update...
File size: 56.3 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer                                           *
4// *                                                                  *
5// * The  Geant4 software  is  copyright of the Copyright Holders  of *
6// * the Geant4 Collaboration.  It is provided  under  the terms  and *
7// * conditions of the Geant4 Software License,  included in the file *
8// * LICENSE and available at  http://cern.ch/geant4/license .  These *
9// * include a list of copyright holders.                             *
10// *                                                                  *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work  make  any representation or  warranty, express or implied, *
14// * regarding  this  software system or assume any liability for its *
15// * use.  Please see the license in the file  LICENSE  and URL above *
16// * for the full disclaimer and the limitation of liability.         *
17// *                                                                  *
18// * This  code  implementation is the result of  the  scientific and *
19// * technical work of the GEANT4 collaboration.                      *
20// * By using,  copying,  modifying or  distributing the software (or *
21// * any work based  on the software)  you  agree  to acknowledge its *
22// * use  in  resulting  scientific  publications,  and indicate your *
23// * acceptance of all terms of the Geant4 Software license.          *
24// ********************************************************************
25//
26//
27// $Id: G4GMocrenFileSceneHandler.cc,v 1.14 2010/05/11 10:44:11 allison Exp $
28// GEANT4 tag $Name:  $
29//
30//
31// Created:  Mar. 31, 2009  Akinori Kimura 
32//           Sep. 22, 2009  Akinori Kimura : modify and fix code to support
33//                                          PrimitiveScorers and clean up
34//
35// GMocrenFile scene handler
36
37
38//----- header files
39#include <fstream>
40#include <cstdlib>
41#include <cstring>
42
43#include "globals.hh"
44#include "G4GMocrenFile.hh"
45#include "G4GMocrenFileSceneHandler.hh"
46#include "G4GMocrenFileViewer.hh"
47#include "G4Point3D.hh"
48#include "G4VisAttributes.hh"
49#include "G4Scene.hh"
50#include "G4Transform3D.hh"
51#include "G4Polyhedron.hh"
52#include "G4Box.hh"
53#include "G4Cons.hh"
54#include "G4Polyline.hh"
55#include "G4Trd.hh"
56#include "G4Tubs.hh"
57#include "G4Trap.hh"
58#include "G4Torus.hh"
59#include "G4Sphere.hh"
60#include "G4Para.hh"
61#include "G4Text.hh"
62#include "G4Circle.hh"
63#include "G4Square.hh"
64#include "G4VPhysicalVolume.hh"
65#include "G4PhysicalVolumeModel.hh"
66#include "G4LogicalVolume.hh"
67#include "G4Material.hh"
68
69#include "G4VPVParameterisation.hh"
70#include "G4VVolumeMaterialScanner.hh"
71#include "G4VisTrajContext.hh"
72#include "G4VisManager.hh"
73#include "G4VTrajectoryModel.hh"
74#include "G4TrajectoryDrawByCharge.hh"
75#include "G4HitsModel.hh"
76#include "G4GMocrenMessenger.hh"
77#include "G4PSHitsModel.hh"
78#include "G4GMocrenIO.hh"
79#include "G4VNestedParameterisation.hh"
80#include "G4GMocrenTouchable.hh"
81#include "G4GMocrenFileCTtoDensityMap.hh"
82#include "G4PhantomParameterisation.hh"
83
84#include "G4ScoringManager.hh"
85#include "G4ScoringBox.hh"
86
87//----- constants
88const char  GDD_FILE_HEADER      [] = "g4_";
89const char  DEFAULT_GDD_FILE_NAME[] = "g4_00.gdd";
90
91const G4int FR_MAX_FILE_NUM = 100 ;
92const G4int MAX_NUM_TRAJECTORIES = 100000;
93
94//-- for a debugging
95const G4bool GFDEBUG = false;
96const G4bool GFDEBUG_TRK = false;//true;
97const G4bool GFDEBUG_HIT = false;//true;
98const G4int GFDEBUG_DET = 0; // 0: false
99
100//////////////////////
101// static variables //
102//////////////////////
103
104//----- static variables
105G4int G4GMocrenFileSceneHandler::kSceneIdCount = 0; 
106
107///////////////////////////
108// Driver-dependent part //
109///////////////////////////
110
111
112//----- G4GMocrenFileSceneHandler, constructor
113G4GMocrenFileSceneHandler::G4GMocrenFileSceneHandler(G4GMocrenFile& system,
114                                                     G4GMocrenMessenger & messenger,
115                                                     const G4String& name)
116  : G4VSceneHandler(system, kSceneIdCount++, name),
117    kSystem(system),
118    kMessenger(messenger),
119    kgMocrenIO(new G4GMocrenIO()),
120    kbSetModalityVoxelSize(false),
121    kbModelingTrajectory(false),
122//    kGddDest(0),
123    kFlagInModeling(false),
124    kFlagSaving_g4_gdd(false),
125    kFlagParameterization(0),
126    kFLagProcessedInteractiveScorer(false) {
127
128  // g4.gdd filename and its directory
129  if(getenv("G4GMocrenFile_DEST_DIR") == NULL) {
130    std::strcpy(kGddDestDir , "");                    // output dir
131    std::strcpy(kGddFileName, DEFAULT_GDD_FILE_NAME); // filename
132  } else {
133    std::strcpy(kGddDestDir , getenv("G4GMocrenFile_DEST_DIR")); // output dir
134    std::strcpy(kGddFileName, DEFAULT_GDD_FILE_NAME); // filename
135  }
136               
137  // maximum number of g4.gdd files in the dest directory
138  kMaxFileNum = FR_MAX_FILE_NUM ; // initialization
139  if ( std::getenv( "G4GMocrenFile_MAX_FILE_NUM" ) != NULL ) { 
140               
141    std::sscanf( getenv("G4GMocrenFile_MAX_FILE_NUM"), "%d", &kMaxFileNum ) ;
142
143  } else {
144    kMaxFileNum = FR_MAX_FILE_NUM ;
145  }
146  if( kMaxFileNum < 1 ) { kMaxFileNum = 1 ; }
147
148  InitializeParameters();
149
150} 
151
152
153//----- G4GMocrenFileSceneHandler, destructor
154G4GMocrenFileSceneHandler::~G4GMocrenFileSceneHandler () 
155{
156  if(GFDEBUG) G4cerr << "***** ~G4GMocrenFileSceneHandler" << G4endl;
157
158  if(kGddDest) {
159    //----- End of modeling
160    // close g4.gdd
161    GFEndModeling();
162  }
163  ClearStore (); // clear current scene
164
165}
166
167//----- initialize all parameters
168void G4GMocrenFileSceneHandler::InitializeParameters() {
169
170  kbSetModalityVoxelSize = false;
171
172  for(G4int i = 0; i < 3; i++) {
173    kModalitySize[i] = 0;
174    kNestedVolumeDimension[i] = 0;
175    kNestedVolumeDirAxis[i] = -1;
176  }
177
178}
179
180//-----
181void G4GMocrenFileSceneHandler::SetGddFileName() 
182{
183  // g4_00.gdd, g4_01.gdd, ..., g4_MAX_FILE_INDEX.gdd
184  const G4int MAX_FILE_INDEX = kMaxFileNum - 1 ;
185
186  // dest directory (null if no environmental variables is set)
187  std::strcpy ( kGddFileName, kGddDestDir) ; 
188
189  // create full path name (default)
190  std::strcat ( kGddFileName, DEFAULT_GDD_FILE_NAME );
191
192  // Automatic updation of file names
193  static G4int currentNumber = 0;
194  for( G4int i = currentNumber ; i < kMaxFileNum ; i++) { 
195
196    // Message in the final execution
197    if( i == MAX_FILE_INDEX ) 
198      {
199        G4cerr << "==========================================="   << G4endl; 
200        G4cerr << "WARNING MESSAGE from GMocrenFile driver:   "   << G4endl;
201        G4cerr << "  This file name is the final one in the   "   << G4endl;
202        G4cerr << "  automatic updation of the output file name." << G4endl; 
203        G4cerr << "  You may overwrite existing files, i.e.   "   << G4endl; 
204        G4cerr << "  g4_XX.gdd."   << G4endl;
205        G4cerr << "==========================================="   << G4endl; 
206      }
207
208    // re-determine file name as G4GMocrenFile_DEST_DIR/g4_XX.gdd
209    if( i >=  0 && i <= 9 ) { 
210      std::sprintf( kGddFileName, "%s%s%s%d.gdd" , kGddDestDir,  GDD_FILE_HEADER, "0", i );
211    } else {
212      std::sprintf( kGddFileName, "%s%s%d.gdd" , kGddDestDir,  GDD_FILE_HEADER, i );
213    }
214
215    // check validity of the file name
216    std::ifstream fin(kGddFileName); 
217    if(GFDEBUG)
218      G4cout << "FILEOPEN: " << i << " : " << kGddFileName << fin.fail() << G4endl;
219    if(!fin) { 
220      // new file       
221      fin.close();
222      currentNumber = i+1;
223      break; 
224    } else { 
225      // already exists (try next)
226      fin.close(); 
227    } 
228
229  } // for
230
231  G4cerr << "======================================================================" << G4endl; 
232  G4cerr << "Output file: " << kGddFileName                          << G4endl; 
233  G4cerr << "Destination directory (current dir if NULL): " << kGddDestDir << G4endl; 
234  G4cerr << "Maximum number of files in the destination directory: " << kMaxFileNum << G4endl; 
235  G4cerr << "Note:" << G4endl; 
236  G4cerr << "  * The maximum number is customizable as:           " << G4endl;
237  G4cerr << "      % setenv  G4GMocrenFile_MAX_FILE_NUM  number " << G4endl;       
238  G4cerr << "  * The destination directory is customizable as:" << G4endl;
239  G4cerr << "      % setenv  G4GMocrenFile_DEST_DIR  dir_name/  " << G4endl;       
240  G4cerr << "     ** Do not forget \"/\" at the end of the dir_name, e.g. \"./tmp/\"." << G4endl;             
241  //G4cerr << "        dir_name, e.g. \"./tmp/\"."                 << G4endl;             
242  G4cerr << G4endl;
243  G4cerr << "Maximum number of trajectories is set to " << MAX_NUM_TRAJECTORIES << "."<< G4endl;
244  G4cerr << "======================================================================" << G4endl; 
245
246} // G4GMocrenFileSceneHandler::SetGddFileName()
247
248
249//-----
250void    G4GMocrenFileSceneHandler::BeginSavingGdd( void )
251{
252  if(GFDEBUG) G4cerr << "***** BeginSavingGdd (called)" << G4endl;
253
254  if( !IsSavingGdd() ) {
255
256    if(GFDEBUG) {
257      G4cerr << "*****                   (started) " ;
258      G4cerr << "(open g4.gdd, ##)"  << G4endl;
259    }
260
261    SetGddFileName() ; // result set to kGddFileName
262    kFlagSaving_g4_gdd = true; 
263
264
265    G4GMocrenFileCTtoDensityMap ctdens;
266    short minmax[2];
267    minmax[0] = ctdens.GetMinCT();
268    minmax[1] = ctdens.GetMaxCT();
269    kgMocrenIO->setModalityImageMinMax(minmax);
270    std::vector<G4float> map;
271    G4float dens;
272    for(G4int i = minmax[0]; i <= minmax[1]; i++) {
273      dens = ctdens.GetDensity(i);
274      map.push_back(dens);
275    }
276    kgMocrenIO->setModalityImageDensityMap(map);
277
278    /*
279    G4String fname = "modality-map.dat";
280    std::ifstream ifile(fname);
281    if(ifile) {
282      short minmax[2];
283      ifile >> minmax[0] >> minmax[1];
284      kgMocrenIO->setModalityImageMinMax(minmax);
285      std::vector<G4float> map;
286      G4float dens;
287      for(G4int i = minmax[0]; i <= minmax[1]; i++) {
288        ifile >> dens;
289        map.push_back(dens);
290      }
291      kgMocrenIO->setModalityImageDensityMap(map);
292     
293    } else {
294      G4cerr << "cann't open the file : " << fname << G4endl;
295    }
296    */
297
298    // mesh size
299    //kMessenger.getNoVoxels(kModalitySize[0], kModalitySize[1], kModalitySize[2]);
300    //kgMocrenIO->setModalityImageSize(kModalitySize);
301   
302    // initializations
303    //kgMocrenIO->clearModalityImage();
304    kgMocrenIO->clearDoseDistAll();
305    kgMocrenIO->clearROIAll();
306    kgMocrenIO->clearTracks();
307    kgMocrenIO->clearDetector();
308    std::vector<Detector>::iterator itr = kDetectors.begin();
309    for(; itr != kDetectors.end(); itr++) {
310      itr->clear();
311    }
312    kDetectors.clear();
313   
314    kNestedHitsList.clear();
315    kNestedVolumeNames.clear();
316     
317  }
318}
319
320void    G4GMocrenFileSceneHandler::EndSavingGdd  ( void ) 
321{
322  if(GFDEBUG) G4cerr << "***** EndSavingGdd (called)" << G4endl;
323
324  if(IsSavingGdd()) {
325    if(GFDEBUG) G4cerr << "*****                 (started) (close "
326                       << kGddFileName << ")" << G4endl;
327
328    if(kGddDest) kGddDest.close();
329    kFlagSaving_g4_gdd = false; 
330
331    std::map<Index3D, G4float>::iterator itr = kNestedModality.begin();
332    G4int xmax=0, ymax=0, zmax=0;
333    for(; itr != kNestedModality.end(); itr++) {
334      if(itr->first.x > xmax) xmax = itr->first.x;
335      if(itr->first.y > ymax) ymax = itr->first.y;
336      if(itr->first.z > zmax) zmax = itr->first.z;
337    }
338    // mesh size
339    kModalitySize[0] = xmax+1;
340    kModalitySize[1] = ymax+1;
341    kModalitySize[2] = zmax+1;
342    kgMocrenIO->setModalityImageSize(kModalitySize);
343    if(GFDEBUG) G4cout << "gMocren-file driver : modality size : "
344                       << kModalitySize[0] << " x "
345                       << kModalitySize[1] << " x "
346                       << kModalitySize[2] << G4endl;
347
348    G4int nxy = kModalitySize[0]*kModalitySize[1];
349    //std::map<G4int, G4float>::iterator itr;
350    for(G4int z = 0; z < kModalitySize[2]; z++) {
351      short * modality = new short[nxy];
352      for(G4int y = 0; y < kModalitySize[1]; y++) {
353        for(G4int x = 0; x < kModalitySize[0]; x++) {
354          //for(G4int x = kModalitySize[0]-1; x >= 0 ; x--) {
355          //G4int ixy = x + (kModalitySize[1]-y-1)*kModalitySize[0];
356
357          G4int ixy = x + y*kModalitySize[0];
358          Index3D idx(x,y,z);
359          itr = kNestedModality.find(idx);
360          if(itr != kNestedModality.end()) {
361
362            modality[ixy] = kgMocrenIO->convertDensityToHU(itr->second);
363          } else {
364            modality[ixy] = -1024;
365          }
366
367        }
368      }
369      kgMocrenIO->setModalityImage(modality);
370    }
371
372    //-- dose
373    size_t nhits = kNestedHitsList.size();
374    if(GFDEBUG) G4cout << "gMocren-file driver : # hits = " << nhits << G4endl;
375
376    std::map<Index3D, G4double>::iterator hitsItr;
377    std::map<G4String, std::map<Index3D, G4double> >::iterator hitsListItr = kNestedHitsList.begin();
378
379    for(G4int n = 0; hitsListItr != kNestedHitsList.end(); hitsListItr++, n++) {
380
381      kgMocrenIO->newDoseDist();
382      kgMocrenIO->setDoseDistName(hitsListItr->first, n);
383      kgMocrenIO->setDoseDistSize(kModalitySize, n);
384
385      G4double minmax[2] = {DBL_MAX, -DBL_MAX};
386      for(G4int z = 0 ; z < kModalitySize[2]; z++) {
387        G4double * values = new G4double[nxy];
388        for(G4int y = 0; y < kModalitySize[1]; y++) {
389          for(G4int x = 0; x < kModalitySize[0]; x++) {
390
391            G4int ixy = x + y*kModalitySize[0];
392            Index3D idx(x,y,z);
393            hitsItr = hitsListItr->second.find(idx);
394            if(hitsItr != hitsListItr->second.end()) {
395
396              values[ixy] = hitsItr->second;
397            } else {
398              values[ixy] = 0.;
399            }
400            if(values[ixy] < minmax[0]) minmax[0] = values[ixy];
401            if(values[ixy] > minmax[1]) minmax[1] = values[ixy];
402          }
403        }
404        kgMocrenIO->setDoseDist(values, n);
405      }
406      kgMocrenIO->setDoseDistMinMax(minmax, n);
407      G4double lower = 0.;
408      if(minmax[0] < 0)  lower = minmax[0];
409      G4double scale = (minmax[1]-lower)/25000.;
410      kgMocrenIO->setDoseDistScale(scale, n);
411      G4String sunit("unit?"); //temporarily
412      kgMocrenIO->setDoseDistUnit(sunit, n);
413    }
414     
415
416    //-- draw axes
417    if(false) {//true,false
418      G4ThreeVector trans;
419      G4RotationMatrix rot;
420      trans = kVolumeTrans3D.getTranslation();
421      rot = kVolumeTrans3D.getRotation().inverse();
422      // x
423      std::vector<G4float *> tracks;
424      unsigned char colors[3];
425      G4float * trk = new G4float[6];
426      tracks.push_back(trk);
427
428      G4ThreeVector orig(0.,0.,0), xa(2000.,0.,0.), ya(0.,2000.,0.), za(0.,0.,2000.);
429      orig -= trans;
430      orig.transform(rot);
431      xa -= trans;
432      xa.transform(rot);
433      ya -= trans;
434      ya.transform(rot);
435      za -= trans;
436      za.transform(rot);
437      for(G4int i = 0; i < 3; i++) trk[i] = orig[i];
438      for(G4int i = 0; i < 3; i++) trk[i+3] = xa[i];
439      colors[0] = 255; colors[1] = 0; colors[2] = 0;
440      kgMocrenIO->addTrack(tracks, colors);
441      // y
442      for(G4int i = 0; i < 3; i++) trk[i+3] = ya[i];
443      colors[0] = 0; colors[1] = 255; colors[2] = 0;
444      kgMocrenIO->addTrack(tracks, colors);
445      // z
446      for(G4int i = 0; i < 3; i++) trk[i+3] = za[i];
447      colors[0] = 0; colors[1] = 0; colors[2] = 255;
448      kgMocrenIO->addTrack(tracks, colors);
449    }
450
451    //-- detector
452    ExtractDetector();
453
454
455    if(GFDEBUG_DET) G4cout << ">>>>>>>>>>>>>>>>>>>>>>   (";
456    std::vector<G4float> transformObjects;
457    for(G4int i = 0; i < 3; i++) {
458      // need to check!!
459      transformObjects.push_back((kVolumeSize[i]/2. - kVoxelDimension[i]/2.));
460      if(GFDEBUG_DET) G4cout << transformObjects[i] << ", ";
461    }
462    if(GFDEBUG_DET) G4cout << ")" << G4endl;
463
464
465    kgMocrenIO->translateTracks(transformObjects);
466    kgMocrenIO->translateDetector(transformObjects);
467
468    // store
469    kgMocrenIO->storeData(kGddFileName);
470  } 
471
472}
473
474
475//-----
476void G4GMocrenFileSceneHandler::GFBeginModeling( void )
477{
478  G4VSceneHandler::BeginModeling();
479
480  if( !GFIsInModeling() ) {
481
482
483      if(GFDEBUG) G4cerr << "***** G4GMocrenFileSceneHandler::GFBeginModeling (called & started)" << G4endl;
484
485      //----- Send saving command and heading comment
486      BeginSavingGdd();
487     
488      kFlagInModeling = true ;
489
490      // These models are entrusted to user commands /vis/scene/add/psHits or hits
491      //scene->AddEndOfEventModel(new G4PSHitsModel());
492      //scene->AddEndOfEventModel(new G4HitsModel());
493      if(GFDEBUG_HIT) {
494        G4Scene * scene = GetScene();
495        std::vector<G4VModel*> vmodel = scene->GetEndOfEventModelList();
496        std::vector<G4VModel*>::iterator itr = vmodel.begin();
497        for(; itr != vmodel.end(); itr++) {
498          G4cout << " IIIIII model name: " << (*itr)->GetGlobalTag() << G4endl;
499        }
500      }
501
502    } // if
503} 
504
505
506//========== AddPrimitive() functions ==========//
507
508//----- Add polyline
509void G4GMocrenFileSceneHandler::AddPrimitive (const G4Polyline& polyline) 
510{
511  if(GFDEBUG) G4cerr << "***** AddPrimitive" << G4endl;
512
513
514  //----- Initialize if necessary
515  GFBeginModeling();
516
517  static G4int numTrajectories = 0;
518  if(numTrajectories >= MAX_NUM_TRAJECTORIES) return;
519
520  // draw trajectories
521  if(kbModelingTrajectory) {
522   
523    G4TrajectoriesModel * pTrModel = dynamic_cast<G4TrajectoriesModel*>(fpModel);
524    if (!pTrModel) { 
525      G4Exception
526        ("G4VSceneHandler::AddCompound(const G4Polyline&): Not a G4TrajectoriesModel.");
527    }
528
529    G4ThreeVector trans;
530    G4RotationMatrix rot;
531    trans = kVolumeTrans3D.getTranslation();
532    rot = kVolumeTrans3D.getRotation().inverse();
533
534    if(GFDEBUG_TRK) G4cout << "   trajectory points : " << G4endl;
535    std::vector<G4float *> trajectory;
536    if(polyline.size() < 2) return;
537    G4Polyline::const_iterator preitr = polyline.begin();
538    G4Polyline::const_iterator postitr = preitr; postitr++;
539    for(; postitr != polyline.end(); preitr++, postitr++) {
540      G4ThreeVector prePts(preitr->x(), preitr->y(), preitr->z());
541      prePts -= trans;
542      prePts.transform(rot);
543      G4ThreeVector postPts(postitr->x(), postitr->y(), postitr->z());
544      postPts -= trans;
545      postPts.transform(rot);
546      G4float * stepPts = new G4float[6];
547      stepPts[0] = prePts.x();
548      stepPts[1] = prePts.y();
549      stepPts[2] = prePts.z();
550      stepPts[3] = postPts.x();
551      stepPts[4] = postPts.y();
552      stepPts[5] = postPts.z();
553      trajectory.push_back(stepPts);
554
555      if(GFDEBUG_TRK) {
556        G4cout << "                       ("
557               << stepPts[0] << ", "
558               << stepPts[1] << ", "
559               << stepPts[2] << ") - (" 
560               << stepPts[3] << ", "
561               << stepPts[4] << ", "
562               << stepPts[5] << ")" << G4endl;
563      }
564    }
565
566    const G4VisAttributes * att = polyline.GetVisAttributes();
567    G4Color color = att->GetColor();
568    unsigned char trkcolor[3];
569    trkcolor[0] = (unsigned char)(color.GetRed()*255);
570    trkcolor[1] = (unsigned char)(color.GetGreen()*255);
571    trkcolor[2] = (unsigned char)(color.GetBlue()*255);
572    if(GFDEBUG_TRK) {
573      G4cout << "              color  : ["
574             << color.GetRed() << ", "
575             << color.GetGreen() << ", "
576             << color.GetBlue() << "]" << G4endl;
577    }
578
579    kgMocrenIO->addTrack(trajectory, trkcolor);
580
581    numTrajectories++;
582  }
583
584} // G4GMocrenFileSceneHandler::AddPrimitive (polyline)
585
586
587//----- Add nurbes
588void G4GMocrenFileSceneHandler::AddPrimitive (const G4NURBS&)
589{
590  //-----
591  if(GFDEBUG) G4cerr << "***** AddPrimitive( G4NURBS )" << G4endl;
592
593  //----- Initialize if necessary
594  GFBeginModeling();
595
596}
597
598
599
600//----- Add text
601void G4GMocrenFileSceneHandler::AddPrimitive ( const G4Text& text )
602{
603  // to avoid a warning in the compile process
604  G4Text dummytext = text;
605
606  //-----
607  if(GFDEBUG) G4cerr << "***** AddPrimitive( G4Text )" << G4endl;
608
609  //----- Initialize IF NECESSARY
610  GFBeginModeling();
611
612} // G4GMocrenFileSceneHandler::AddPrimitive ( text )
613
614
615//----- Add circle
616void G4GMocrenFileSceneHandler::AddPrimitive ( const G4Circle& mark_circle )
617{
618  // to avoid a warning in the compile process
619  G4Circle dummycircle = mark_circle;
620
621  //-----
622  if(GFDEBUG) G4cerr << "***** AddPrimitive( G4Circle )" << G4endl;
623
624  //----- Initialize IF NECESSARY
625  GFBeginModeling();
626
627
628} // G4GMocrenFileSceneHandler::AddPrimitive ( mark_circle )
629
630
631//----- Add square
632void G4GMocrenFileSceneHandler::AddPrimitive (const G4Square& mark_square )
633{
634  // to avoid a warning in the compile process
635  G4Square dummysquare = mark_square;
636
637  //-----
638  if(GFDEBUG) G4cerr << "***** AddPrimitive( G4Square )" << G4endl;
639
640  //----- Initialize if necessary
641  GFBeginModeling();
642
643} // G4GMocrenFileSceneHandler::AddPrimitive ( mark_square )
644
645
646//----- Add polyhedron
647void G4GMocrenFileSceneHandler::AddPrimitive ( const G4Polyhedron& polyhedron ) 
648{
649  //-----
650  if(GFDEBUG) G4cerr << "***** AddPrimitive( G4Polyhedron )" << G4endl;
651
652
653  if (polyhedron.GetNoFacets() == 0) return;
654
655  //----- Initialize if necessary
656  GFBeginModeling();
657
658  //---------- (3) Facet block
659  for (G4int f = polyhedron.GetNoFacets(); f; f--){
660    G4int notLastEdge;
661    G4int index = -1; // initialization
662    G4int edgeFlag = 1;
663    G4int preedgeFlag = 1;
664    G4int work[4], i = 0;
665    do {
666      preedgeFlag = edgeFlag;
667      notLastEdge = polyhedron.GetNextVertexIndex(index, edgeFlag);
668      work[i++] = index;
669    }while (notLastEdge);
670    switch (i){
671    case 3:
672      //SendStrInt3(FR_FACET, work[0], work[1], work[2] );
673      break;
674    case 4:
675      //SendStrInt4(FR_FACET, work[0], work[1], work[2], work[3] );
676      break;
677    default:
678      G4cerr <<
679        "ERROR G4GMocrenFileSceneHandler::AddPrimitive(G4Polyhedron)" << G4endl;
680      G4PhysicalVolumeModel* pPVModel =
681        dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
682      if (pPVModel) G4cerr <<
683                      "Volume " << pPVModel->GetCurrentPV()->GetName() <<
684                      ", Solid " << pPVModel->GetCurrentLV()->GetSolid()->GetName() <<
685                      " (" << pPVModel->GetCurrentLV()->GetSolid()->GetEntityType();
686      G4cerr <<
687        "\nG4Polyhedron facet with " << i << " edges" << G4endl;       
688    }
689  }
690
691} // G4GMocrenFileSceneHandler::AddPrimitive (polyhedron)
692
693
694//-----
695void G4GMocrenFileSceneHandler::GFEndModeling ()
696{
697  G4VSceneHandler::EndModeling();
698
699  //-----               
700  if(GFDEBUG) G4cerr << "***** GFEndModeling (called)" << G4endl;
701
702  if( GFIsInModeling() ) {
703
704    if(GFDEBUG) {
705      G4cerr << "***** GFEndModeling (started) " ; 
706      G4cerr << "(/EndModeling, /DrawAll, /CloseDevice)" << G4endl;
707    }
708
709    //----- End saving data to g4.gdd
710    EndSavingGdd() ;
711
712    //------ Reset flag
713    kFlagInModeling = false ;
714
715  }
716
717}
718
719
720//-----
721void G4GMocrenFileSceneHandler::BeginPrimitives (const G4Transform3D& objectTransformation)
722{
723  if(GFDEBUG) G4cerr << "***** BeginPrimitives " << G4endl;
724
725  GFBeginModeling();
726
727  G4VSceneHandler::BeginPrimitives (objectTransformation);
728  fpObjectTransformation = &objectTransformation;
729
730}
731
732
733//-----
734void G4GMocrenFileSceneHandler::EndPrimitives ()
735{
736  if(GFDEBUG) G4cerr << "***** EndPrimitives " << G4endl;
737
738  G4VSceneHandler::EndPrimitives ();
739}
740
741
742//========== AddSolid() functions ==========//
743
744//----- Add box
745void G4GMocrenFileSceneHandler::AddSolid( const G4Box& box )
746{
747  if(GFDEBUG) G4cerr << "***** AddSolid ( box )" << G4endl;
748
749  if(GFDEBUG_DET > 0)
750    G4cout << "G4GMocrenFileSceneHandler::AddSolid(const G4Box&)  : "
751           << box.GetName() << G4endl;
752
753  //----- skip drawing invisible primitive
754  if( !IsVisible() ) { return ; }
755
756  //----- Initialize if necessary
757  GFBeginModeling();
758
759
760  //--
761  if(GFDEBUG_DET > 1) {
762    G4cout << "-------" << G4endl;
763    G4cout << "    " << box.GetName() << G4endl;
764    G4Polyhedron * poly = box.CreatePolyhedron();
765    poly->Transform(*fpObjectTransformation);
766    //G4int nv = poly->GetNoVertices();
767    G4Point3D v1, v2;
768    G4int next;
769    //while(1) { // next flag isn't functional.
770    for(G4int i = 0; i < 12; i++) { // # of edges is 12.
771      poly->GetNextEdge(v1, v2, next);
772      if(next == 0) break;
773      G4cout << "    (" << v1.x() << ", "
774             << v1.y() << ", "
775             << v1.z() << ") - ("
776             << v2.x() << ", "
777             << v2.y() << ", "
778             << v2.z() << ") [" << next << "]"
779             << G4endl;
780    }
781    delete poly;
782  }
783
784
785  // the volume name set by /vis/gMocren/setVolumeName
786  G4String volName = kMessenger.getVolumeName();
787
788  if(kFlagParameterization != 2) {
789    G4ScoringManager * pScrMan = G4ScoringManager::GetScoringManager();
790    if(pScrMan) {
791      G4ScoringBox * pScBox = dynamic_cast<G4ScoringBox*>(pScrMan->FindMesh(volName));
792      G4bool bMesh = false;
793      if(pScBox != NULL) bMesh = true;
794      if(bMesh) kFlagParameterization = 2;
795      if(GFDEBUG_DET > 0) G4cout << "   G4ScoringManager::FindMesh() : "
796                                 << volName << " - " << bMesh << G4endl;
797    }
798  }
799
800  const G4VModel* pv_model  = GetModel();
801  if (!pv_model) { return ; } 
802  G4PhysicalVolumeModel* pPVModel =
803    dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
804  if (!pPVModel) { return ; }
805
806
807  //-- debug information
808  if(GFDEBUG_DET > 0) {
809    G4Material * mat = pPVModel->GetCurrentMaterial();
810    G4String name = mat->GetName();
811    G4double dens = mat->GetDensity()/(g/cm3);
812    G4int copyNo = pPVModel->GetCurrentPV()->GetCopyNo();
813    G4int depth = pPVModel->GetCurrentDepth();
814    G4cout << "    copy no.: " << copyNo << G4endl;
815    G4cout << "    depth   : " << depth << G4endl;
816    G4cout << "    density : " << dens << " [g/cm3]" << G4endl;
817    G4cout << "    location: " << pPVModel->GetCurrentPV()->GetObjectTranslation() << G4endl;
818    G4cout << "    Multiplicity        : " << pPVModel->GetCurrentPV()->GetMultiplicity() << G4endl;
819    G4cout << "    Is replicated?      : " << pPVModel->GetCurrentPV()->IsReplicated() << G4endl;
820    G4cout << "    Is parameterised?   : " << pPVModel->GetCurrentPV()->IsParameterised() << G4endl;
821    G4cout << "    top phys. vol. name : " << pPVModel->GetTopPhysicalVolume()->GetName() << G4endl;
822  }
823
824  //-- check the parameterised volume
825  if(box.GetName() == volName) {
826   
827    kVolumeTrans3D = *fpObjectTransformation;
828    // coordination system correction for gMocren
829    G4ThreeVector raxis(1., 0., 0.), dummy(0.,0.,0.); 
830    G4RotationMatrix rot(raxis, pi*rad);
831    G4Transform3D trot(rot, dummy);
832    if(GFDEBUG_DET) {
833      G4ThreeVector trans = kVolumeTrans3D.getTranslation();
834      G4RotationMatrix rot = kVolumeTrans3D.getRotation().inverse();
835      G4cout << "kVolumeTrans3D: " << trans << G4endl << rot << G4endl;
836    }
837    kVolumeTrans3D = kVolumeTrans3D*trot;
838    if(GFDEBUG_DET) G4cout << " Parameterised volume : " << box.GetName() << G4endl;
839
840
841
842    //
843    G4VPhysicalVolume * pv[3] = {0,0,0};
844    pv[0] = pPVModel->GetCurrentPV()->GetLogicalVolume()->GetDaughter(0);
845    if(!pv[0]) {
846      G4Exception("Error[gMocrenFileSceneHandler]: Unexpected volume.");
847    }
848    G4int dirAxis[3] = {-1,-1,-1};
849    G4int nDaughters[3] = {0,0,0};
850
851    EAxis axis; G4int nReplicas; G4double width; G4double offset; G4bool consuming; 
852    pv[0]->GetReplicationData(axis, nReplicas, width, offset, consuming);
853    nDaughters[0] = nReplicas;
854    switch(axis) {
855    case kXAxis: dirAxis[0] = 0; break;
856    case kYAxis: dirAxis[0] = 1; break;
857    case kZAxis: dirAxis[0] = 2; break;
858    default: G4Exception("Error.");
859    }
860    kNestedVolumeNames.push_back(pv[0]->GetName());
861    if(GFDEBUG_DET) 
862      G4cout << "        daughter name :  " << pv[0]->GetName()
863             << "   # : " << nDaughters[0] << G4endl;
864
865    //
866    if(GFDEBUG_DET) {
867      if(pv[0]->GetLogicalVolume()->GetNoDaughters()) {
868        G4cout << "# of daughters : " 
869               << pv[0]->GetLogicalVolume()->GetNoDaughters() << G4endl;
870      } else {
871        //G4Exception("Error: G4GMocrenFileSceneHandler: 0000010");
872      }
873    }
874
875    // check whether nested or regular parameterization
876    if(GFDEBUG_DET) G4cout << "# of daughters : " 
877                           << pv[0]->GetLogicalVolume()->GetNoDaughters() << G4endl;
878    if(pv[0]->GetLogicalVolume()->GetNoDaughters() == 0) {
879      kFlagParameterization = 1;
880      //G4Exception("Error: G4GMocrenFileSceneHandler: 0000020");
881    }
882   
883    if(kFlagParameterization == 0) {
884
885      pv[1] = pv[0]->GetLogicalVolume()->GetDaughter(0);
886      if(pv[1]) {
887        pv[1]->GetReplicationData(axis, nReplicas, width, offset, consuming);
888        nDaughters[1] = nReplicas;
889        switch(axis) {
890        case kXAxis: dirAxis[1] = 0; break;
891        case kYAxis: dirAxis[1] = 1; break;
892        case kZAxis: dirAxis[1] = 2; break;
893        default: G4Exception("Error.");
894        }
895        kNestedVolumeNames.push_back(pv[1]->GetName());
896        if(GFDEBUG_DET) 
897          G4cout << "        sub-daughter name :  " << pv[1]->GetName()
898                 << "   # : " << nDaughters[1]<< G4endl;
899
900        //
901        pv[2] = pv[1]->GetLogicalVolume()->GetDaughter(0);
902        if(pv[2]) {
903          nDaughters[2] = pv[2]->GetMultiplicity();
904          kNestedVolumeNames.push_back(pv[2]->GetName());
905          if(GFDEBUG_DET) 
906            G4cout << "        sub-sub-daughter name :  " << pv[2]->GetName()
907                   << "   # : " << nDaughters[2] << G4endl;
908
909          if(nDaughters[2] > 1) {
910            G4VNestedParameterisation * nestPara
911              = dynamic_cast<G4VNestedParameterisation*>(pv[2]->GetParameterisation());
912            if(!nestPara) G4Exception("Error[gMocrenFileSceneHandler]: None nested parameterisation");
913            nestPara->ComputeTransformation(0, pv[2]);
914            G4ThreeVector trans0 = pv[2]->GetObjectTranslation();
915            nestPara->ComputeTransformation(1, pv[2]);
916            G4ThreeVector trans1 = pv[2]->GetObjectTranslation();
917            G4ThreeVector diff(trans0 - trans1);
918            if(GFDEBUG_DET) 
919              G4cout << trans0 << " - " << trans1 << " - " << diff << G4endl;
920
921            if(diff.x() != 0.) dirAxis[2] = 0;
922            else if(diff.y() != 0.) dirAxis[2] = 1;
923            else if(diff.z() != 0.) dirAxis[2] = 2;
924            else G4Exception("Error[gMocrenFileSceneHandler]: Unexpected nested parameterisation");
925          }
926        }
927      }
928     
929      for(G4int i = 0; i < 3; i++) {
930        kNestedVolumeDimension[i] = nDaughters[i];
931        //kNestedVolumeDimension[i] = nDaughters[dirAxis[i]];
932        kNestedVolumeDirAxis[i] = dirAxis[i];
933      }
934      //G4cout << "@@@@@@@@@ "
935      //       << dirAxis[0] << ", " << dirAxis[1] << ", " << dirAxis[2] << G4endl;
936
937      // get densities
938      G4VNestedParameterisation * nestPara
939        = dynamic_cast<G4VNestedParameterisation*>(pv[2]->GetParameterisation());
940      if(nestPara) {
941        G4double prexyz[3] = {0.,0.,0.}, xyz[3] = {0.,0.,0.};
942        for(G4int n0 = 0; n0 < nDaughters[0]; n0++) {
943          for(G4int n1 = 0; n1 < nDaughters[1]; n1++) {
944            for(G4int n2 = 0; n2 < nDaughters[2]; n2++) {
945                 
946              G4GMocrenTouchable * touch = new G4GMocrenTouchable(n1, n0);
947              if(GFDEBUG_DET) 
948                G4cout << "   retrieve volume : copy # : " << n0
949                       << ", " << n1 << ", " << n2 << G4endl;
950              G4Material * mat = nestPara->ComputeMaterial(pv[2], n2, touch);
951              delete touch;
952              G4double dens = mat->GetDensity()/(g/cm3);
953 
954              if(GFDEBUG_DET) 
955                G4cout << "           density :" << dens << " [g/cm3]" << G4endl;
956
957              G4Box tbox(box);
958              nestPara->ComputeDimensions(tbox, n2, pv[2]);
959              xyz[0] = tbox.GetXHalfLength()/mm;
960              xyz[1] = tbox.GetYHalfLength()/mm;
961              xyz[2] = tbox.GetZHalfLength()/mm;
962              if(n0 != 0 || n1 != 0 || n2 != 0) {
963                for(G4int i = 0; i < 3; i++) {
964                  if(xyz[i] != prexyz[i]) G4Exception("Error[gMocrenFileSceneHandler]: Unsupported parameterisation.");
965                }
966              }
967              if(GFDEBUG_DET) 
968                G4cout << "              size : " << tbox.GetXHalfLength()/mm << " x "
969                       << tbox.GetYHalfLength()/mm << " x "
970                       << tbox.GetZHalfLength()/mm << " [mm3]" << G4endl;
971
972              G4int idx[3];
973              idx[dirAxis[0]] = n0;
974              idx[dirAxis[1]] = n1;
975              idx[dirAxis[2]] = n2;
976              Index3D i3d(idx[0],idx[1],idx[2]);
977              kNestedModality[i3d] = dens;
978              if(GFDEBUG_DET) 
979                G4cout << " index: " << idx[0] << ", " << idx[1] << ", " << idx[2]
980                       << "  density: " << dens << G4endl;
981
982              for(G4int i = 0; i < 3; i++) prexyz[i] = xyz[i];
983            }
984          }
985        } 
986
987        kVolumeSize.set(box.GetXHalfLength()*2/mm,
988                        box.GetYHalfLength()*2/mm,
989                        box.GetZHalfLength()*2/mm);
990        // mesh size
991        if(!kbSetModalityVoxelSize) {
992          G4float spacing[3] = {static_cast<G4float>(2*xyz[0]),
993                                static_cast<G4float>(2*xyz[1]),
994                                static_cast<G4float>(2*xyz[2])};
995          kgMocrenIO->setVoxelSpacing(spacing);
996          kVoxelDimension.set(spacing[0], spacing[1], spacing[2]);
997          kbSetModalityVoxelSize = true;
998        }
999
1000      } else {
1001        if(GFDEBUG_DET) 
1002          G4cout << pv[2]->GetName() << G4endl;
1003        G4Exception("Error[gMocrenFileSceneHandler]: none nested parameterization");
1004      }
1005
1006
1007
1008      //-- debug
1009      if(GFDEBUG_DET > 1) {
1010        if(pPVModel->GetCurrentPV()->IsParameterised()) {
1011          G4VPVParameterisation * para = pPVModel->GetCurrentPV()->GetParameterisation();
1012          G4cout << " Is nested parameterisation? : " << para->IsNested() << G4endl;
1013
1014
1015          G4int npvp = pPVModel->GetDrawnPVPath().size();
1016          G4cout << "     physical volume node id : "
1017                 << "size: " << npvp << ", PV name: ";
1018          for(G4int i = 0; i < npvp; i++) {
1019            G4cout << pPVModel->GetDrawnPVPath()[i].GetPhysicalVolume()->GetName()
1020                   << " [param:"
1021                   << pPVModel->GetDrawnPVPath()[i].GetPhysicalVolume()->IsParameterised()
1022                   << ",rep:"
1023                   << pPVModel->GetDrawnPVPath()[i].GetPhysicalVolume()->IsReplicated();
1024            if(pPVModel->GetDrawnPVPath()[i].GetPhysicalVolume()->GetParameterisation()) {
1025              G4cout << ",nest:"
1026                     << pPVModel->GetDrawnPVPath()[i].GetPhysicalVolume()->GetParameterisation()->IsNested();
1027            }
1028            G4cout << ",copyno:"
1029                   << pPVModel->GetDrawnPVPath()[i].GetPhysicalVolume()->GetCopyNo();
1030            G4cout << "] - ";
1031          }
1032          G4cout << G4endl;
1033
1034
1035          EAxis axis; G4int nReplicas; G4double width; G4double offset; G4bool consuming; 
1036          pPVModel->GetCurrentPV()->GetReplicationData(axis, nReplicas, width, offset, consuming);
1037          G4cout << "     # replicas : " << nReplicas << G4endl;
1038          G4double pareDims[3] = {0.,0.,0.};
1039          G4Box * pbox = dynamic_cast<G4Box *>(pPVModel->GetDrawnPVPath()[npvp-2].GetPhysicalVolume()->GetLogicalVolume()->GetSolid());
1040          if(pbox) {
1041            pareDims[0] = 2.*pbox->GetXHalfLength()*mm;
1042            pareDims[1] = 2.*pbox->GetYHalfLength()*mm;
1043            pareDims[2] = 2.*pbox->GetZHalfLength()*mm;
1044            G4cout << "     mother size ["
1045                   << pPVModel->GetDrawnPVPath()[npvp-2].GetPhysicalVolume()->GetName()
1046                   << "] : "
1047                   << pareDims[0] << " x "
1048                   << pareDims[1] << " x "
1049                   << pareDims[2] << " [mm3]"
1050                   << G4endl;
1051          }
1052          G4double paraDims[3];
1053          G4Box * boxP = dynamic_cast<G4Box *>(pPVModel->GetDrawnPVPath()[npvp-1].GetPhysicalVolume()->GetLogicalVolume()->GetSolid());
1054          if(boxP) {
1055            paraDims[0] = 2.*boxP->GetXHalfLength()*mm;
1056            paraDims[1] = 2.*boxP->GetYHalfLength()*mm;
1057            paraDims[2] = 2.*boxP->GetZHalfLength()*mm;
1058            G4cout << "     parameterised volume? ["
1059                   << pPVModel->GetDrawnPVPath()[npvp-1].GetPhysicalVolume()->GetName()
1060                   << "] : "
1061                   << paraDims[0] << " x "
1062                   << paraDims[1] << " x "
1063                   << paraDims[2] << " [mm3]  : "
1064                   << G4int(pareDims[0]/paraDims[0]) << " x "
1065                   << G4int(pareDims[1]/paraDims[1]) << " x "
1066                   << G4int(pareDims[2]/paraDims[2]) << G4endl;
1067          } else {
1068            G4cout << pPVModel->GetDrawnPVPath()[npvp-2].GetPhysicalVolume()->GetName()
1069                   << " isn't a G4Box." << G4endl;
1070          }
1071        }
1072      }
1073
1074
1075    } else if(kFlagParameterization == 1) { // G4PhantomParameterisation based geom. construnction
1076
1077      // get the dimension of the parameterized patient geometry
1078      G4PhantomParameterisation * phantomPara
1079        = dynamic_cast<G4PhantomParameterisation*>(pv[0]->GetParameterisation());
1080      if(!phantomPara) {
1081        G4Exception("Error: G4GMocrenFileSceneHandler: no G4PhantomParameterisation");
1082      } else {
1083        ;
1084      }
1085
1086      kNestedVolumeDimension[0] = phantomPara->GetNoVoxelX();
1087      kNestedVolumeDimension[1] = phantomPara->GetNoVoxelY();
1088      kNestedVolumeDimension[2] = phantomPara->GetNoVoxelZ();
1089      kNestedVolumeDirAxis[0] = 0;
1090      kNestedVolumeDirAxis[1] = 1;
1091      kNestedVolumeDirAxis[2] = 2;
1092
1093      // get densities of the parameterized patient geometry
1094      G4int nX = kNestedVolumeDimension[0];
1095      G4int nXY = kNestedVolumeDimension[0]*kNestedVolumeDimension[1];
1096
1097      for(G4int n0 = 0; n0 < kNestedVolumeDimension[0]; n0++) {
1098        for(G4int n1 = 0; n1 < kNestedVolumeDimension[1]; n1++) {
1099          for(G4int n2 = 0; n2 < kNestedVolumeDimension[2]; n2++) {
1100
1101            G4int repNo = n0 + n1*nX + n2*nXY;
1102            G4Material * mat = phantomPara->ComputeMaterial(repNo, pv[0]);
1103            G4double dens = mat->GetDensity()/(g/cm3);
1104
1105
1106            G4int idx[3];
1107            idx[kNestedVolumeDirAxis[0]] = n0;
1108            idx[kNestedVolumeDirAxis[1]] = n1;
1109            idx[kNestedVolumeDirAxis[2]] = n2;
1110            Index3D i3d(idx[0],idx[1],idx[2]);
1111            kNestedModality[i3d] = dens;
1112
1113            if(GFDEBUG_DET) 
1114              G4cout << " index: " << idx[0] << ", " << idx[1] << ", " << idx[2]
1115                     << "  density: " << dens << G4endl;
1116
1117          }
1118        }
1119      }
1120
1121      kVolumeSize.set(box.GetXHalfLength()*2/mm,
1122                      box.GetYHalfLength()*2/mm,
1123                      box.GetZHalfLength()*2/mm);
1124
1125      // mesh size
1126      if(!kbSetModalityVoxelSize) {
1127        G4float spacing[3] = {static_cast<G4float>(2*phantomPara->GetVoxelHalfX()),
1128                              static_cast<G4float>(2*phantomPara->GetVoxelHalfY()),
1129                              static_cast<G4float>(2*phantomPara->GetVoxelHalfZ())};
1130        kgMocrenIO->setVoxelSpacing(spacing);
1131        kVoxelDimension.set(spacing[0], spacing[1], spacing[2]);
1132        kbSetModalityVoxelSize = true;
1133      }
1134    }
1135
1136  } // if(box.GetName() == volName)
1137
1138
1139  // processing geometry construction based on the interactive PS
1140  if(!kFLagProcessedInteractiveScorer) {
1141
1142
1143    // get the dimension of the geometry defined in G4VScoringMesh
1144    G4ScoringManager * pScrMan = G4ScoringManager::GetScoringManager();
1145    if(!pScrMan) return;
1146    G4ScoringBox * scoringBox
1147      = dynamic_cast<G4ScoringBox*>(pScrMan->FindMesh(volName));
1148    if(scoringBox == NULL) return;
1149
1150
1151    G4int nVoxels[3];
1152    scoringBox->GetNumberOfSegments(nVoxels);
1153    // this order depends on the G4ScoringBox
1154    kNestedVolumeDimension[0] = nVoxels[2];
1155    kNestedVolumeDimension[1] = nVoxels[1];
1156    kNestedVolumeDimension[2] = nVoxels[0];
1157    kNestedVolumeDirAxis[0] = 2;
1158    kNestedVolumeDirAxis[1] = 1;
1159    kNestedVolumeDirAxis[2] = 0;
1160
1161    // get densities of the parameterized patient geometry
1162    for(G4int n0 = 0; n0 < kNestedVolumeDimension[0]; n0++) {
1163      for(G4int n1 = 0; n1 < kNestedVolumeDimension[1]; n1++) {
1164        for(G4int n2 = 0; n2 < kNestedVolumeDimension[2]; n2++) {
1165
1166          G4double dens = 0.*(g/cm3);
1167
1168          G4int idx[3];
1169          idx[kNestedVolumeDirAxis[0]] = n0;
1170          idx[kNestedVolumeDirAxis[1]] = n1;
1171          idx[kNestedVolumeDirAxis[2]] = n2;
1172          Index3D i3d(idx[0],idx[1],idx[2]);
1173          kNestedModality[i3d] = dens;
1174
1175        }
1176      }
1177    }
1178
1179    G4ThreeVector boxSize = scoringBox->GetSize();
1180    if(GFDEBUG_DET > 1) {
1181      G4cout << "Interactive Scorer : size - "
1182             << boxSize.x()/cm << " x " 
1183             << boxSize.y()/cm << " x " 
1184             << boxSize.z()/cm << " [cm3]" << G4endl;
1185      G4cout << "Interactive Scorer : # voxels - "
1186             << nVoxels[0] << " x " 
1187             << nVoxels[1] << " x " 
1188             << nVoxels[2] << G4endl;
1189    }
1190    kVolumeSize.set(boxSize.x()*2,
1191                    boxSize.y()*2,
1192                    boxSize.z()*2);
1193
1194    // mesh size
1195    if(!kbSetModalityVoxelSize) {
1196      G4float spacing[3] = {static_cast<G4float>(boxSize.x()*2/nVoxels[0]),
1197                            static_cast<G4float>(boxSize.y()*2/nVoxels[1]),
1198                            static_cast<G4float>(boxSize.z()*2/nVoxels[2])};
1199
1200      kgMocrenIO->setVoxelSpacing(spacing);
1201      kVoxelDimension.set(spacing[0], spacing[1], spacing[2]);
1202      kbSetModalityVoxelSize = true;
1203
1204    }
1205
1206
1207    kVolumeTrans3D = *fpObjectTransformation;
1208
1209    // translation for the scoring mesh
1210    G4ThreeVector sbth = scoringBox->GetTranslation();
1211    G4Translate3D sbtranslate(sbth);
1212    kVolumeTrans3D = kVolumeTrans3D*sbtranslate;
1213
1214    // rotation matrix for the scoring mesh
1215    G4RotationMatrix sbrm;
1216    sbrm = scoringBox->GetRotationMatrix();
1217    if(!sbrm.isIdentity()) {
1218      G4ThreeVector sbdummy(0.,0.,0.); 
1219      G4Transform3D sbrotate(sbrm.inverse(), sbdummy);
1220      kVolumeTrans3D = kVolumeTrans3D*sbrotate;
1221    }
1222
1223
1224    // coordination system correction for gMocren
1225    G4ThreeVector raxisY(0., 1., 0.), dummyY(0.,0.,0.); 
1226    G4RotationMatrix rotY(raxisY, pi*rad);
1227    G4Transform3D trotY(rotY, dummyY);
1228    G4ThreeVector raxisZ(0., 0., 1.), dummyZ(0.,0.,0.); 
1229    G4RotationMatrix rotZ(raxisZ, pi*rad);
1230    G4Transform3D trotZ(rotZ, dummyZ);
1231
1232    kVolumeTrans3D = kVolumeTrans3D*trotY*trotZ;
1233
1234
1235    //
1236    kFLagProcessedInteractiveScorer = true;
1237  } 
1238
1239
1240
1241  //-- add detectors
1242  G4bool bAddDet = true;
1243  if(!kMessenger.getDrawVolumeGrid()) {
1244
1245    if(kFlagParameterization == 0) { // nested parameterisation
1246
1247      if(volName == box.GetName()) {
1248        bAddDet = false;
1249      }
1250
1251      std::vector<G4String>::iterator itr = kNestedVolumeNames.begin();
1252      for(; itr != kNestedVolumeNames.end(); itr++) {
1253        if(*itr == box.GetName())  {
1254          bAddDet = false;
1255          break;
1256        }
1257      }
1258    } else if(kFlagParameterization == 1) { // phantom paramemterisation
1259
1260      if(volName != box.GetName()) {
1261        bAddDet = false;
1262      }
1263
1264    } else if(kFlagParameterization == 2) { // interactive primitive scorer
1265      ;
1266    }
1267
1268  }
1269  if(bAddDet) AddDetector(box);
1270
1271
1272} // void G4GMocrenFileSceneHandler::AddSolid( const G4Box& box )
1273
1274
1275//----- Add tubes
1276void 
1277G4GMocrenFileSceneHandler::AddSolid( const G4Tubs& tubes )
1278{
1279  if(GFDEBUG) G4cerr << "***** AddSolid ( tubes )" << G4endl;
1280
1281  //----- skip drawing invisible primitive
1282  if( !IsVisible() ) { return ; }
1283
1284  //----- Initialize if necessary
1285  GFBeginModeling();
1286
1287  //
1288  AddDetector(tubes);
1289
1290
1291  // for a debug
1292  if(GFDEBUG_DET > 0) {
1293    G4cout << "-------" << G4endl;
1294    G4cout << "    " << tubes.GetName() << G4endl;
1295    G4Polyhedron * poly = tubes.CreatePolyhedron();
1296    G4int nv = poly->GetNoVertices();
1297    for(G4int i = 0; i < nv; i++) {
1298      G4cout << "    (" << poly->GetVertex(i).x() << ", "
1299             << poly->GetVertex(i).y() << ", "
1300             << poly->GetVertex(i).z() << ")" << G4endl;
1301    }
1302    delete poly;
1303  }
1304
1305  const G4VModel* pv_model  = GetModel();
1306  if (!pv_model) { return ; } 
1307  G4PhysicalVolumeModel* pPVModel =
1308    dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
1309  if (!pPVModel) { return ; }
1310  G4Material * mat = pPVModel->GetCurrentMaterial();
1311  G4String name = mat->GetName();
1312
1313} // void G4GMocrenFileSceneHandler::AddSolid( const G4Tubs& )
1314
1315
1316
1317//----- Add cons
1318void 
1319G4GMocrenFileSceneHandler::AddSolid( const G4Cons& cons )
1320{
1321  if(GFDEBUG) G4cerr << "***** AddSolid ( cons )" << G4endl;
1322
1323  //----- skip drawing invisible primitive
1324  if( !IsVisible() ) { return ; }
1325
1326  //----- Initialize if necessary
1327  GFBeginModeling();
1328
1329  //
1330  AddDetector(cons);
1331
1332}// G4GMocrenFileSceneHandler::AddSolid( cons )
1333
1334
1335//----- Add trd
1336void G4GMocrenFileSceneHandler::AddSolid ( const G4Trd& trd )
1337{
1338  if(GFDEBUG) G4cerr << "***** AddSolid ( trd )" << G4endl;
1339
1340
1341  //----- skip drawing invisible primitive
1342  if( !IsVisible() ) { return ; }
1343
1344  //----- Initialize if necessary
1345  GFBeginModeling();
1346
1347  //
1348  AddDetector(trd);
1349
1350} // G4GMocrenFileSceneHandler::AddSolid ( trd )
1351
1352
1353//----- Add sphere
1354void G4GMocrenFileSceneHandler::AddSolid ( const G4Sphere& sphere )
1355{
1356  if(GFDEBUG) G4cerr << "***** AddSolid ( sphere )" << G4endl;
1357
1358  //----- skip drawing invisible primitive
1359  if( !IsVisible() ) { return ; }
1360
1361  //----- Initialize if necessary
1362  GFBeginModeling();
1363
1364  //
1365  AddDetector(sphere);
1366
1367} // G4GMocrenFileSceneHandler::AddSolid ( sphere )
1368
1369
1370//----- Add para
1371void G4GMocrenFileSceneHandler::AddSolid (const G4Para& para)
1372{
1373  if(GFDEBUG) G4cerr << "***** AddSolid ( para )" << G4endl;
1374
1375  //----- skip drawing invisible primitive
1376  if( !IsVisible() ) { return ; }
1377
1378  //----- Initialize if necessary
1379  GFBeginModeling();
1380
1381  //
1382  AddDetector(para);
1383
1384} // G4GMocrenFileSceneHandler::AddSolid ( para )
1385
1386
1387//----- Add trap
1388void G4GMocrenFileSceneHandler::AddSolid (const G4Trap& trap)
1389{
1390  if(GFDEBUG) G4cerr << "***** AddSolid ( trap )" << G4endl;
1391
1392  //----- skip drawing invisible primitive
1393  if( !IsVisible() ) { return ; }
1394
1395  //----- Initialize if necessary
1396  GFBeginModeling();
1397
1398  //
1399  AddDetector(trap);
1400
1401} // G4GMocrenFileSceneHandler::AddSolid (const G4Trap& trap)
1402
1403
1404//----- Add torus
1405void 
1406G4GMocrenFileSceneHandler::AddSolid( const G4Torus& torus )
1407{
1408  if(GFDEBUG) G4cerr << "***** AddSolid ( torus )" << G4endl;
1409
1410  //----- skip drawing invisible primitive
1411  if( !IsVisible() ) { return ; }
1412
1413  //----- Initialize if necessary
1414  GFBeginModeling();
1415
1416  //
1417  AddDetector(torus);
1418
1419} // void G4GMocrenFileSceneHandler::AddSolid( const G4Torus& )
1420
1421
1422
1423//----- Add a shape which is not treated above
1424void G4GMocrenFileSceneHandler::AddSolid ( const G4VSolid& solid  )
1425{
1426  //----- skip drawing invisible primitive
1427  if( !IsVisible() ) { return ; }
1428
1429  //----- Initialize if necessary
1430  GFBeginModeling();
1431
1432  //
1433  AddDetector(solid);
1434
1435  //----- Send a primitive
1436  G4VSceneHandler::AddSolid( solid ) ; 
1437
1438} //G4GMocrenFileSceneHandler::AddSolid ( const G4VSolid& )
1439
1440#include "G4TrajectoriesModel.hh"
1441#include "G4VTrajectory.hh"
1442#include "G4VTrajectoryPoint.hh"
1443
1444//----- Add a trajectory
1445void G4GMocrenFileSceneHandler::AddCompound(const G4VTrajectory & traj) {
1446
1447  kbModelingTrajectory = true;
1448
1449  G4VSceneHandler::AddCompound(traj);
1450
1451  if(GFDEBUG_TRK) {
1452    std::cout << " ::AddCompound(const G4VTrajectory&) >>>>>>>>> " << std::endl;
1453    G4TrajectoriesModel * pTrModel = dynamic_cast<G4TrajectoriesModel*>(fpModel);
1454    if (!pTrModel) { 
1455      G4Exception
1456        ("G4VSceneHandler::AddCompound(const G4VTrajectory&): Not a G4TrajectoriesModel.");
1457    } else {
1458      traj.DrawTrajectory(pTrModel->GetDrawingMode());
1459
1460      const G4VTrajectory * trj = pTrModel->GetCurrentTrajectory();
1461      G4cout << "------ track" << G4endl;
1462      G4cout << "    name:     " << trj->GetParticleName() << G4endl;
1463      G4cout << "    id:       " << trj->GetTrackID() << G4endl;
1464      G4cout << "    charge:   " << trj->GetCharge() << G4endl;
1465      G4cout << "    momentum: " << trj->GetInitialMomentum() << G4endl;
1466     
1467      G4int nPnt = trj->GetPointEntries();
1468      G4cout << "    point:    ";
1469      for(G4int i = 0; i < nPnt; i++) {
1470        G4cout << trj->GetPoint(i)->GetPosition() << ", ";
1471      }
1472      G4cout << G4endl;
1473    }
1474    G4cout << G4endl;
1475  }
1476
1477  kbModelingTrajectory = false;
1478}
1479
1480#include <vector>
1481#include "G4VHit.hh"
1482#include "G4AttValue.hh"
1483//----- Add a hit
1484void G4GMocrenFileSceneHandler::AddCompound( const G4VHit & hit) {
1485  if(GFDEBUG_HIT) G4cout << " ::AddCompound(const G4VHit&) >>>>>>>>> " << G4endl;
1486
1487  G4VSceneHandler::AddCompound(hit);
1488
1489  /*
1490    const std::map<G4String, G4AttDef> * map = hit.GetAttDefs();
1491    if(!map) return;
1492    std::map<G4String, G4AttDef>::const_iterator itr = map->begin();
1493    for(; itr != map->end(); itr++) {
1494    G4cout << itr->first << " : " << itr->second.GetName()
1495    << " , " << itr->second.GetDesc() << G4endl;
1496    }
1497  */
1498
1499  std::vector<G4String> hitNames = kMessenger.getHitNames();
1500  if(GFDEBUG_HIT) {
1501    std::vector<G4String>::iterator itr = hitNames.begin();
1502    for(; itr != hitNames.end(); itr++) 
1503      G4cout << "  hit name : " << *itr << G4endl;
1504  }
1505 
1506  std::vector<G4AttValue> * attval = hit.CreateAttValues();
1507  if(!attval) {G4cout << "0 empty " << (unsigned long)attval << G4endl;}
1508  else {
1509
1510    G4bool bid[3] = {false, false, false};
1511    Index3D id;
1512
1513    std::vector<G4AttValue>::iterator itr;
1514    // First, get IDs
1515    for(itr = attval->begin(); itr != attval->end(); itr++) {
1516      std::string stmp = itr->GetValue();
1517      std::istringstream sval(stmp.c_str());
1518
1519      if(itr->GetName() == G4String("XID")) {
1520        sval >> id.x;
1521        bid[0] = true;
1522        continue;
1523      }
1524      if(itr->GetName() == G4String("YID")) {
1525        sval >> id.y;
1526        bid[1] = true;
1527        continue;
1528      }
1529      if(itr->GetName() == G4String("ZID")) {
1530        sval >> id.z;
1531        bid[2] = true;
1532        continue;
1533      }
1534    }
1535
1536    G4int nhitname = (G4int)hitNames.size();
1537
1538    if(bid[0] && bid[1] && bid[2]) {
1539
1540      if(GFDEBUG_HIT)
1541        G4cout << " Hit : index(" << id.x << ", " << id.y << ", "
1542               << id.z << ")" << G4endl;
1543
1544      // Get attributes
1545      for(itr = attval->begin(); itr != attval->end(); itr++) {
1546        for(G4int i = 0; i < nhitname; i++) {
1547          if(itr->GetName() == hitNames[i]) {
1548
1549            std::string stmp = itr->GetValue();
1550            std::istringstream sval(stmp.c_str());
1551            G4double value;
1552            G4String unit;
1553            sval >> value >> unit;
1554
1555            std::map<G4String, std::map<Index3D, G4double> >::iterator kNestedHitsListItr;
1556            kNestedHitsListItr = kNestedHitsList.find(hitNames[i]);
1557            if(kNestedHitsListItr != kNestedHitsList.end()) {
1558              //fTempNestedHits = &kNestedHitsListItr->second;
1559              //(*fTempNestedHits)[id] = value;
1560              kNestedHitsListItr->second[id] = value;
1561            } else {
1562              std::map<Index3D, G4double> hits;
1563              hits.insert(std::map<Index3D, G4double>::value_type(id, value));
1564              kNestedHitsList[hitNames[i]] = hits;
1565            }
1566
1567           
1568            if(GFDEBUG_HIT)
1569              G4cout << "     : " << hitNames[i] << " -> " << value
1570                     << " [" << unit << "]" << G4endl;
1571          }
1572        }
1573      }
1574    } else {
1575      G4Exception("Error in G4GMocrenFileSceneHandler::AddCompound(const G4VHit &)");
1576    }
1577
1578    delete attval;
1579  }
1580
1581}
1582
1583void G4GMocrenFileSceneHandler::AddCompound(const G4THitsMap<G4double> & hits) {
1584  if(GFDEBUG_HIT)
1585    G4cout << " ::AddCompound(const std::map<G4int, G4double*> &) >>>>>>>>> " << G4endl;
1586
1587
1588  std::vector<G4String> hitScorerNames = kMessenger.getHitScorerNames();
1589  G4int nhitname = (G4int)hitScorerNames.size();
1590  G4String scorername = static_cast<G4VHitsCollection>(hits).GetName();
1591
1592  //-- --//
1593  /*
1594  std::vector<G4String> hitScorerNames = kMessenger.getHitScorerNames();
1595  if(GFDEBUG_HIT) {
1596    std::vector<G4String>::iterator itr = hitScorerNames.begin();
1597    for(; itr != hitScorerNames.end(); itr++)
1598      G4cout << "  PS name : " << *itr << G4endl;
1599  }
1600  */
1601 
1602
1603  //for(G4int i = 0; i < nhitname; i++) {       // this selection trusts
1604    //if(scorername == hitScorerNames[i]) {   // thea command /vis/scene/add/psHits hit_name.
1605
1606      G4int idx[3];
1607      std::map<G4int, G4double*> * map = hits.GetMap();
1608      std::map<G4int, G4double*>::const_iterator itr = map->begin();
1609      for(; itr != map->end(); itr++) {
1610        GetNestedVolumeIndex(itr->first, idx);
1611        Index3D id(idx[0], idx[1], idx[2]);
1612       
1613        std::map<G4String, std::map<Index3D, G4double> >::iterator nestedHitsListItr;
1614        nestedHitsListItr = kNestedHitsList.find(scorername);
1615        if(nestedHitsListItr != kNestedHitsList.end()) {
1616          nestedHitsListItr->second[id] = *itr->second;
1617        } else {
1618          std::map<Index3D, G4double> hit;
1619          hit.insert(std::map<Index3D, G4double>::value_type(id, *itr->second));
1620          kNestedHitsList[scorername] = hit;
1621        }
1622      }
1623 
1624      //break;
1625    //}
1626  //}
1627
1628  if(GFDEBUG_HIT) {
1629    G4String meshname = static_cast<G4VHitsCollection>(hits).GetSDname();
1630    G4cout << "       >>>>> " << meshname << " : " << scorername  << G4endl;
1631
1632    for(G4int i = 0; i < nhitname; i++)
1633      if(scorername == hitScorerNames[i]) 
1634        G4cout << "       !!!! Hit scorer !!!! " << scorername << G4endl;
1635
1636    G4cout << " dimension: "
1637           << kNestedVolumeDimension[0] << " x "
1638           << kNestedVolumeDimension[1] << " x "
1639           << kNestedVolumeDimension[2] << G4endl;
1640
1641    G4int id[3];
1642    std::map<G4int, G4double*> * map = hits.GetMap();
1643    std::map<G4int, G4double*>::const_iterator itr = map->begin();
1644    for(; itr != map->end(); itr++) {
1645      GetNestedVolumeIndex(itr->first, id);
1646      G4cout << "[" << itr->first << "] "
1647             << "("<< id[0] << "," << id[1] << "," << id[2] << ")"
1648             << *itr->second << ", ";
1649    }
1650    G4cout << G4endl;
1651  }
1652
1653}
1654
1655//-----
1656G4bool G4GMocrenFileSceneHandler::IsVisible()
1657{
1658  //-----
1659  G4bool  visibility  = true ;
1660
1661  const G4VisAttributes* pVisAttribs =
1662    fpViewer->GetApplicableVisAttributes( fpVisAttribs );
1663
1664  if(pVisAttribs) {
1665    visibility = pVisAttribs->IsVisible();
1666  } 
1667
1668  return visibility ;
1669
1670} // G4GMocrenFileSceneHandler::IsVisible()
1671
1672
1673//-----
1674void G4GMocrenFileSceneHandler::ClearTransientStore() 
1675{
1676  G4VSceneHandler::ClearTransientStore ();
1677  // This is typically called after an update and before drawing hits
1678  // of the next event.  To simulate the clearing of "transients"
1679  // (hits, etc.) the detector is redrawn...
1680  if (fpViewer) {
1681    fpViewer -> SetView ();
1682    fpViewer -> ClearView ();
1683    fpViewer -> DrawView ();
1684  }
1685}
1686
1687//-----
1688void G4GMocrenFileSceneHandler::AddDetector(const G4VSolid & solid) {
1689
1690  Detector detector;
1691
1692  // detector name
1693  detector.name = solid.GetName();
1694  if(GFDEBUG_DET > 1)
1695    G4cout << "0 Detector name : " << detector.name << G4endl;
1696 
1697  const G4VModel* pv_model  = GetModel();
1698  if (!pv_model) { return ; } 
1699  G4PhysicalVolumeModel* pPVModel =
1700    dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
1701  if (!pPVModel) { return ; }
1702
1703  // edge points of the detector
1704  std::vector<G4float *> dedges;
1705  G4Polyhedron * poly = solid.CreatePolyhedron();
1706  detector.polyhedron = poly;
1707  detector.transform3D = *fpObjectTransformation;
1708
1709  // retrieve color
1710  unsigned char uccolor[3] = {30, 30, 30};
1711  if(pPVModel->GetCurrentLV()->GetVisAttributes()) {
1712    G4Color color = pPVModel->GetCurrentLV()->GetVisAttributes()->GetColor();
1713    uccolor[0] = (unsigned char)(color.GetRed()*255);
1714    uccolor[1] = (unsigned char)(color.GetGreen()*255);
1715    uccolor[2] = (unsigned char)(color.GetBlue()*255);
1716    //if(uccolor[0] < 2 && uccolor[1] < 2 && uccolor[2] < 2)
1717    //uccolor[0] = uccolor[1] = uccolor[2] = 30; // dark grey
1718  }
1719  for(G4int i = 0; i < 3; i++) detector.color[i] = uccolor[i];
1720  //
1721  kDetectors.push_back(detector);
1722
1723  if(GFDEBUG_DET > 1) {
1724    G4cout << "0     color:   (" << (G4int)uccolor[0] << ", "
1725           << (G4int)uccolor[1] << ", " << (G4int)uccolor[2] << ")"
1726           << G4endl;
1727  }
1728
1729}
1730
1731//-----
1732void G4GMocrenFileSceneHandler::ExtractDetector() {
1733
1734  std::vector<Detector>::iterator itr = kDetectors.begin();
1735
1736  for(; itr != kDetectors.end(); itr++) {
1737
1738    // detector name
1739    G4String detname = itr->name;
1740    if(GFDEBUG_DET > 1)
1741      G4cout << "Detector name : " << detname << G4endl;
1742
1743    // edge points of the detector
1744    std::vector<G4float *> dedges;
1745    G4Polyhedron * poly = itr->polyhedron;
1746    poly->Transform(itr->transform3D);
1747    G4Transform3D invVolTrans = kVolumeTrans3D.inverse();
1748    poly->Transform(invVolTrans);
1749
1750    G4Point3D v1, v2;
1751    G4bool bnext = true;
1752    G4int next;
1753    G4int nedges = 0;
1754    //
1755    while(bnext) {
1756      if(!(poly->GetNextEdge(v1, v2, next))) bnext =false;
1757      G4float * edge = new G4float[6];
1758      edge[0] = v1.x()/mm;
1759      edge[1] = v1.y()/mm;
1760      edge[2] = v1.z()/mm;
1761      edge[3] = v2.x()/mm;
1762      edge[4] = v2.y()/mm;
1763      edge[5] = v2.z()/mm;
1764      dedges.push_back(edge);
1765      nedges++;
1766    }
1767    //delete poly;
1768    // detector color
1769    unsigned char uccolor[3] = {itr->color[0],
1770                                itr->color[1],
1771                                itr->color[2]};
1772    //
1773    kgMocrenIO->addDetector(detname, dedges, uccolor);
1774    for(G4int i = 0; i < nedges; i++) { // # of edges is 12.
1775      delete [] dedges[i];
1776    }
1777    dedges.clear(); 
1778
1779    if(GFDEBUG_DET > 1) {
1780      G4cout << "    color:   (" << (G4int)uccolor[0] << ", "
1781             << (G4int)uccolor[1] << ", " << (G4int)uccolor[2] << ")"
1782             << G4endl;
1783    }
1784  }
1785}
1786
1787void G4GMocrenFileSceneHandler::GetNestedVolumeIndex(G4int _idx, G4int _idx3d[3]) {
1788  if(kNestedVolumeDimension[0] == 0 ||
1789     kNestedVolumeDimension[1] == 0 ||
1790     kNestedVolumeDimension[2] == 0) {
1791    for(G4int i = 0; i < 3; i++) _idx3d[i] = 0;
1792    return;
1793  }
1794
1795
1796  if(kFlagParameterization == 0) {
1797
1798    G4int plane = kNestedVolumeDimension[2]*kNestedVolumeDimension[1];
1799    G4int line = kNestedVolumeDimension[2];
1800 
1801  /*
1802  G4int idx3d[3];
1803  idx3d[0] = _idx/plane;
1804  idx3d[1] = (_idx%plane)/line;
1805  idx3d[2] = (_idx%plane)%line;
1806  _idx3d[0] = idx3d[kNestedVolumeDirAxis[0]];
1807  _idx3d[1] = idx3d[kNestedVolumeDirAxis[1]];
1808  _idx3d[2] = idx3d[kNestedVolumeDirAxis[2]];
1809  */
1810
1811    _idx3d[kNestedVolumeDirAxis[0]] = _idx/plane;
1812    _idx3d[kNestedVolumeDirAxis[1]] = (_idx%plane)/line;
1813    _idx3d[kNestedVolumeDirAxis[2]] = (_idx%plane)%line;
1814
1815
1816
1817  /*
1818
1819  G4cout << "G4GMocrenFileSceneHandler::GetNestedVolumeIndex : " << G4endl;
1820  G4cout << "(depi, depj, depk) : "
1821         << kNestedVolumeDirAxis[0] << ", "
1822         << kNestedVolumeDirAxis[1] << ", "
1823         << kNestedVolumeDirAxis[2] << G4endl;
1824  G4cout << "(ni, nj, nk) :"
1825         << kNestedVolumeDimension[0] << ", "
1826         << kNestedVolumeDimension[1] << ", "
1827         << kNestedVolumeDimension[2] << " - " << G4endl;
1828
1829  G4cout << " _idx = " << _idx << "  :  plane = "
1830         << plane << " ,   line = " << line << G4endl;
1831  G4cout << "(idx,idy,idz) + " << _idx3d[0] << ", "
1832         << _idx3d[1] << ", " << _idx3d[2] << " + " << G4endl;
1833
1834  */
1835
1836
1837
1838  } else {
1839   
1840    G4int plane = kNestedVolumeDimension[0]*kNestedVolumeDimension[1];
1841    G4int line = kNestedVolumeDimension[0];
1842    _idx3d[kNestedVolumeDirAxis[2]] = _idx/plane;
1843    _idx3d[kNestedVolumeDirAxis[1]] = (_idx%plane)/line;
1844    _idx3d[kNestedVolumeDirAxis[0]] = (_idx%plane)%line;
1845
1846  }
1847
1848}
1849
1850
1851//-- --//
1852G4GMocrenFileSceneHandler::Detector::Detector()
1853  : polyhedron(0) {
1854  color[0] = color[1] = color[2] = 255;
1855}
1856G4GMocrenFileSceneHandler::Detector::~Detector() {
1857  if(!polyhedron) delete polyhedron;
1858}
1859void G4GMocrenFileSceneHandler::Detector::clear() {
1860  name.clear();
1861  if(!polyhedron) delete polyhedron;
1862  color[0] = color[1] = color[2] = 255;
1863  transform3D = G4Transform3D::Identity;
1864}
1865
1866//-- --//
1867G4GMocrenFileSceneHandler::Index3D::Index3D()
1868  : x(0), y(0), z(0) {
1869  ;
1870}
1871
1872G4GMocrenFileSceneHandler::Index3D::Index3D(const Index3D & _index3D) 
1873  : x(_index3D.x), y(_index3D.y), z(_index3D.z) {
1874  //: x(_index3D.X()),
1875  //y(_index3D.Y()),
1876  //z(_index3D.Z()) {
1877  //  : x(static_cast<Index3D>(_index3D).x),
1878  //    y(static_cast<Index3D>(_index3D).y),
1879  //    z(static_cast<Index3D>(_index3D).z) {
1880  ;
1881}
1882
1883G4GMocrenFileSceneHandler::Index3D::Index3D(G4int _x, G4int _y, G4int _z) 
1884  : x(_x), y(_y), z(_z) {
1885  ;
1886}
1887G4bool G4GMocrenFileSceneHandler::Index3D::operator < (const Index3D & _right) const {
1888  if(z < static_cast<Index3D>(_right).z) {
1889     return true;
1890  } else if(z == _right.z) {
1891    if(y < static_cast<Index3D>(_right).y) return true;
1892    else if(y == _right.y) 
1893      if(x < static_cast<Index3D>(_right).x) return true;
1894  } 
1895  return false;
1896}
1897G4bool G4GMocrenFileSceneHandler::Index3D::operator == (const Index3D & _right) const {
1898  if(z == _right.z && y == _right.y && x == _right.x) return true;
1899  return false;
1900}
Note: See TracBrowser for help on using the repository browser.