source: trunk/source/visualization/RayTracer/src/G4TheRayTracer.cc @ 1347

Last change on this file since 1347 was 1347, checked in by garnier, 13 years ago

geant4 tag 9.4

File size: 14.0 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: G4TheRayTracer.cc,v 1.3 2006/06/29 21:24:25 gunter Exp $
28// GEANT4 tag $Name: geant4-09-04-ref-00 $
29//
30//
31//
32
33
34#include "G4TheRayTracer.hh"
35#include "G4EventManager.hh"
36#include "G4RTMessenger.hh"
37#include "G4RayShooter.hh"
38#include "G4VFigureFileMaker.hh"
39#include "G4RTTrackingAction.hh"
40#include "G4RTSteppingAction.hh"
41#include "G4RayTrajectory.hh"
42#include "G4RayTrajectoryPoint.hh"
43#include "G4RTJpegMaker.hh"
44#include "G4RTSimpleScanner.hh"
45#include "G4GeometryManager.hh"
46#include "G4SDManager.hh"
47#include "G4StateManager.hh"
48#include "G4Event.hh"
49#include "G4TrajectoryContainer.hh"
50#include "G4Colour.hh"
51#include "G4VisAttributes.hh"
52#include "G4UImanager.hh"
53#include "G4TransportationManager.hh"
54#include "G4RegionStore.hh"
55#include "G4ProductionCutsTable.hh"
56
57G4TheRayTracer::G4TheRayTracer(G4VFigureFileMaker* figMaker,
58                               G4VRTScanner* scanner)
59{
60  theFigMaker = figMaker;
61  if(!theFigMaker) theFigMaker = new G4RTJpegMaker;
62  theScanner = scanner;
63  if(!theScanner) theScanner = new G4RTSimpleScanner;
64  theRayShooter = new G4RayShooter();
65  theRayTracerEventAction = 0;
66  theRayTracerStackingAction = 0;
67  theRayTracerTrackingAction = 0;
68  theRayTracerSteppingAction = 0;
69  theMessenger = G4RTMessenger::GetInstance(this,theRayTracerSteppingAction);
70  theEventManager = G4EventManager::GetEventManager();
71
72  nColumn = 640;
73  nRow = 640;
74
75  eyePosition = G4ThreeVector(1.*m,1.*m,1.*m);
76  targetPosition = G4ThreeVector(0.,0.,0.);
77  lightDirection = G4ThreeVector(-0.1,-0.2,-0.3).unit();
78  viewSpan = 5.0*deg;
79  headAngle = 270.*deg; 
80  attenuationLength = 1.0*m;
81
82  distortionOn = false;
83  antialiasingOn = false;
84
85  backgroundColour = G4Colour(1.,1.,1.);
86}
87
88G4TheRayTracer::~G4TheRayTracer()
89{
90  delete theRayShooter;
91  if(theRayTracerTrackingAction) delete theRayTracerTrackingAction;
92  if(theRayTracerSteppingAction) delete theRayTracerSteppingAction;
93  delete theMessenger;
94  delete theScanner;
95  delete theFigMaker;
96}
97
98void G4TheRayTracer::Trace(G4String fileName)
99{
100  G4StateManager* theStateMan = G4StateManager::GetStateManager();
101  G4ApplicationState currentState = theStateMan->GetCurrentState();
102  if(currentState!=G4State_Idle)
103  {
104    G4cerr << "Illegal application state - Trace() ignored." << G4endl;
105    return;
106  }
107
108  if(!theFigMaker)
109  {
110    G4cerr << "Figure file maker class is not specified - Trace() ignored." << G4endl;
111    return;
112  }
113
114  G4UImanager* UI = G4UImanager::GetUIpointer();
115  G4int storeTrajectory = UI->GetCurrentIntValue("/tracking/storeTrajectory");
116  if(storeTrajectory==0) UI->ApplyCommand("/tracking/storeTrajectory 1");
117
118
119  G4ThreeVector tmpVec = targetPosition - eyePosition;
120  eyeDirection = tmpVec.unit();
121  colorR = new unsigned char[nColumn*nRow];
122  colorG = new unsigned char[nColumn*nRow];
123  colorB = new unsigned char[nColumn*nRow];
124
125  StoreUserActions();
126  G4bool succeeded = CreateBitMap();
127  if(succeeded)
128  { CreateFigureFile(fileName); }
129  else
130  { G4cerr << "Could not create figure file" << G4endl;
131    G4cerr << "You might set the eye position outside of the world volume" << G4endl; }
132  RestoreUserActions();
133
134  if(storeTrajectory==0) UI->ApplyCommand("/tracking/storeTrajectory 0");
135
136  delete [] colorR;
137  delete [] colorG;
138  delete [] colorB;
139}
140
141void G4TheRayTracer::StoreUserActions()
142{ 
143  theUserEventAction = theEventManager->GetUserEventAction();
144  theUserStackingAction = theEventManager->GetUserStackingAction();
145  theUserTrackingAction = theEventManager->GetUserTrackingAction();
146  theUserSteppingAction = theEventManager->GetUserSteppingAction();
147
148  if(!theRayTracerTrackingAction) theRayTracerTrackingAction = new G4RTTrackingAction();
149  if(!theRayTracerSteppingAction) theRayTracerSteppingAction = new G4RTSteppingAction();
150
151  theEventManager->SetUserAction(theRayTracerEventAction);
152  theEventManager->SetUserAction(theRayTracerStackingAction);
153  theEventManager->SetUserAction(theRayTracerTrackingAction);
154  theEventManager->SetUserAction(theRayTracerSteppingAction);
155
156  G4SDManager* theSDMan = G4SDManager::GetSDMpointerIfExist();
157  if(theSDMan)
158  { theSDMan->Activate("/",false); }
159
160  G4GeometryManager* theGeomMan = G4GeometryManager::GetInstance();
161  theGeomMan->OpenGeometry();
162  theGeomMan->CloseGeometry(true);
163}
164
165void G4TheRayTracer::RestoreUserActions()
166{
167  theEventManager->SetUserAction(theUserEventAction);
168  theEventManager->SetUserAction(theUserStackingAction);
169  theEventManager->SetUserAction(theUserTrackingAction);
170  theEventManager->SetUserAction(theUserSteppingAction);
171
172  G4SDManager* theSDMan = G4SDManager::GetSDMpointerIfExist();
173  if(theSDMan)
174  { theSDMan->Activate("/",true); }
175}
176
177#include "G4ProcessManager.hh"
178#include "G4ProcessVector.hh"
179#include "G4Geantino.hh"
180
181G4bool G4TheRayTracer::CreateBitMap()
182{
183  G4int iEvent = 0;
184  G4double stepAngle = viewSpan/100.;
185  G4double viewSpanX = stepAngle*nColumn;
186  G4double viewSpanY = stepAngle*nRow;
187  G4bool succeeded;
188
189// Confirm process(es) of Geantino is initialized
190  G4VPhysicalVolume* pWorld =
191        G4TransportationManager::GetTransportationManager()->
192        GetNavigatorForTracking()->GetWorldVolume();
193  G4RegionStore::GetInstance()->UpdateMaterialList(pWorld);
194  G4ProductionCutsTable::GetProductionCutsTable()->UpdateCoupleTable(pWorld);
195  G4ProcessVector* pVector
196    = G4Geantino::GeantinoDefinition()->GetProcessManager()->GetProcessList();
197  for (G4int j=0; j < pVector->size(); ++j) {
198      (*pVector)[j]->BuildPhysicsTable(*(G4Geantino::GeantinoDefinition()));
199  }
200
201// Close geometry and set the application state
202  G4GeometryManager* geomManager = G4GeometryManager::GetInstance();
203  geomManager->OpenGeometry();
204  geomManager->CloseGeometry(1,0);
205 
206  G4ThreeVector center(0,0,0);
207  G4Navigator* navigator =
208      G4TransportationManager::GetTransportationManager()->GetNavigatorForTracking();
209  navigator->LocateGlobalPointAndSetup(center,0,false);
210
211  G4StateManager* theStateMan = G4StateManager::GetStateManager();
212  theStateMan->SetNewState(G4State_GeomClosed); 
213
214// Event loop
215  theScanner->Initialize(nRow,nColumn);
216  G4int iRow, iColumn;
217  while (theScanner->Coords(iRow,iColumn)) {
218      G4int iCoord = iRow * nColumn + iColumn;
219      G4double dRow = 0, dColumn = 0;  // Antialiasing increments.
220      G4Event* anEvent = new G4Event(iEvent++);
221      G4double angleX = -(viewSpanX/2. - (iColumn+dColumn)*stepAngle);
222      G4double angleY = viewSpanY/2. - (iRow+dRow)*stepAngle;
223      G4ThreeVector rayDirection;
224      if(distortionOn)
225      {
226        rayDirection = G4ThreeVector(std::tan(angleX)/std::cos(angleY),-std::tan(angleY)/std::cos(angleX),1.0);
227      }
228      else
229      {
230        rayDirection = G4ThreeVector(std::tan(angleX),-std::tan(angleY),1.0);
231      }
232      rayDirection.rotateZ(headAngle);
233      rayDirection.rotateUz(eyeDirection);
234      G4ThreeVector rayPosition(eyePosition);
235      G4bool interceptable = true;
236      // Check if rayPosition is in the world.
237      EInside whereisit =
238        pWorld->GetLogicalVolume()->GetSolid()->Inside(rayPosition);
239      if (whereisit != kInside) {
240        // It's outside the world, so move it inside.
241        G4double outsideDistance =
242          pWorld->GetLogicalVolume()->GetSolid()->
243          DistanceToIn(rayPosition,rayDirection); 
244        if (outsideDistance != kInfinity) {
245          // Borrowing from geometry, where 1e-8 < epsilon < 1e-3, in
246          // absolute/internal length units, is used for ensuring good
247          // behaviour, choose to add 0.001 to ensure rayPosition is
248          // definitely inside the world volume (JA 16/9/2005)...
249          rayPosition = rayPosition+(outsideDistance+0.001)*rayDirection;
250        }
251        else {
252          interceptable = false;
253        }
254      }
255      if (interceptable) {
256        theRayShooter->Shoot(anEvent,rayPosition,rayDirection);
257        theEventManager->ProcessOneEvent(anEvent);
258        succeeded = GenerateColour(anEvent);
259        colorR[iCoord] = (unsigned char)(int(255*rayColour.GetRed()));
260        colorG[iCoord] = (unsigned char)(int(255*rayColour.GetGreen()));
261        colorB[iCoord] = (unsigned char)(int(255*rayColour.GetBlue()));
262      } else {  // Ray does not intercept world at all.
263        // Store background colour...
264        colorR[iCoord] = (unsigned char)(int(255*backgroundColour.GetRed()));
265        colorG[iCoord] = (unsigned char)(int(255*backgroundColour.GetGreen()));
266        colorB[iCoord] = (unsigned char)(int(255*backgroundColour.GetBlue()));
267        succeeded = true;
268      }
269
270      theScanner->Draw(colorR[iCoord],colorG[iCoord],colorB[iCoord]);
271
272      delete anEvent;
273      if(!succeeded) return false;
274  }
275
276  theStateMan->SetNewState(G4State_Idle); 
277  return true;
278}
279
280void G4TheRayTracer::CreateFigureFile(G4String fileName)
281{
282  //G4cout << nColumn << " " << nRow << G4endl;
283  theFigMaker->CreateFigureFile(fileName,nColumn,nRow,colorR,colorG,colorB);
284}
285
286G4bool G4TheRayTracer::GenerateColour(G4Event* anEvent)
287{
288  G4TrajectoryContainer * trajectoryContainer = anEvent->GetTrajectoryContainer();
289 
290  G4RayTrajectory* trajectory = (G4RayTrajectory*)( (*trajectoryContainer)[0] );
291  if(!trajectory) return false;
292
293  G4int nPoint = trajectory->GetPointEntries();
294  if(nPoint==0) return false;
295
296  G4Colour initialColour(backgroundColour);
297  if( trajectory->GetPointC(nPoint-1)->GetPostStepAtt() )
298  { initialColour = GetSurfaceColour(trajectory->GetPointC(nPoint-1)); }
299  rayColour = Attenuate(trajectory->GetPointC(nPoint-1),initialColour);
300
301  for(int i=nPoint-2;i>=0;i--)
302  {
303    G4Colour surfaceColour = GetSurfaceColour(trajectory->GetPointC(i));
304    G4double weight = 1.0 - surfaceColour.GetAlpha();
305    G4Colour mixedColour = GetMixedColour(rayColour,surfaceColour,weight);
306    rayColour = Attenuate(trajectory->GetPointC(i),mixedColour);
307  }
308   
309  return true;
310}
311
312G4Colour G4TheRayTracer::GetMixedColour(G4Colour surfCol,G4Colour transCol,G4double weight)
313{
314  G4double r = weight*surfCol.GetRed() + (1.-weight)*transCol.GetRed();
315  G4double g = weight*surfCol.GetGreen() + (1.-weight)*transCol.GetGreen();
316  G4double b = weight*surfCol.GetBlue() + (1.-weight)*transCol.GetBlue();
317  G4double a = weight*surfCol.GetAlpha() + (1.-weight)*transCol.GetAlpha();
318  return G4Colour(r,g,b,a);
319}
320
321G4Colour G4TheRayTracer::GetSurfaceColour(G4RayTrajectoryPoint* point)
322{
323  const G4VisAttributes* preAtt = point->GetPreStepAtt();
324  const G4VisAttributes* postAtt = point->GetPostStepAtt();
325
326  G4bool preVis = ValidColour(preAtt);
327  G4bool postVis = ValidColour(postAtt);
328
329  G4Colour transparent(1.,1.,1.,0.);
330
331  if(!preVis&&!postVis) return transparent;
332
333  G4ThreeVector normal = point->GetSurfaceNormal();
334
335  G4Colour preCol(1.,1.,1.);
336  G4Colour postCol(1.,1.,1.);
337
338  if(preVis)
339  {
340    G4double brill = (1.0-(-lightDirection).dot(normal))/2.0;
341    G4double r = preAtt->GetColour().GetRed();
342    G4double g = preAtt->GetColour().GetGreen();
343    G4double b = preAtt->GetColour().GetBlue();
344    preCol = G4Colour(r*brill,g*brill,b*brill,preAtt->GetColour().GetAlpha());
345  }
346  else
347  { preCol = transparent; }
348
349  if(postVis)
350  {
351    G4double brill = (1.0-(-lightDirection).dot(-normal))/2.0;
352    G4double r = postAtt->GetColour().GetRed();
353    G4double g = postAtt->GetColour().GetGreen();
354    G4double b = postAtt->GetColour().GetBlue();
355    postCol = G4Colour(r*brill,g*brill,b*brill,postAtt->GetColour().GetAlpha());
356  }
357  else
358  { postCol = transparent; }
359   
360  if(!preVis) return postCol;
361  if(!postVis) return preCol;
362
363  G4double weight = 0.5;
364  return GetMixedColour(preCol,postCol,weight);
365}
366
367G4Colour G4TheRayTracer::Attenuate(G4RayTrajectoryPoint* point, G4Colour sourceCol)
368{
369  const G4VisAttributes* preAtt = point->GetPreStepAtt();
370
371  G4bool visible = ValidColour(preAtt);
372  if(!visible) return sourceCol;
373
374  G4Colour objCol = preAtt->GetColour();
375  G4double stepRed = objCol.GetRed();
376  G4double stepGreen = objCol.GetGreen();
377  G4double stepBlue = objCol.GetBlue();
378  G4double stepAlpha = objCol.GetAlpha();
379  G4double stepLength = point->GetStepLength();
380
381  G4double attenuationFuctor;
382  if(stepAlpha > 0.9999999){ stepAlpha = 0.9999999; } // patch to the next line
383    attenuationFuctor = -stepAlpha/(1.0-stepAlpha)*stepLength/attenuationLength;
384 
385  G4double KtRed = std::exp((1.0-stepRed)*attenuationFuctor);
386  G4double KtGreen = std::exp((1.0-stepGreen)*attenuationFuctor);
387  G4double KtBlue = std::exp((1.0-stepBlue)*attenuationFuctor);
388  if(KtRed>1.0){KtRed=1.0;}
389  if(KtGreen>1.0){KtGreen=1.0;}
390  if(KtBlue>1.0){KtBlue=1.0;}
391  return G4Colour(sourceCol.GetRed()*KtRed,
392    sourceCol.GetGreen()*KtGreen,sourceCol.GetBlue()*KtBlue);
393}
394
395G4bool G4TheRayTracer::ValidColour(const G4VisAttributes* visAtt)
396{
397  G4bool val = true;
398  if(!visAtt)
399  { val = false; }
400  else if(!(visAtt->IsVisible()))
401  { val = false; }
402  else if(visAtt->IsForceDrawingStyle()
403    &&(visAtt->GetForcedDrawingStyle()==G4VisAttributes::wireframe))
404  { val = false; }
405  return val;
406}
407
Note: See TracBrowser for help on using the repository browser.