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