ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VoxelLimits.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4VoxelLimits.cc
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 // Class G4VoxelLimits implementation
27 //
28 // 13.07.95, P.Kent - Initial version
29 // --------------------------------------------------------------------
30 
31 #include "G4VoxelLimits.hh"
32 
33 #include "G4ios.hh"
34 
36 //
37 // Empty constructor and destructor
38 //
40 {
41 }
42 
44 {
45 }
46 
48 //
49 // Further restrict limits
50 // No checks for illegal restrictions
51 //
52 void G4VoxelLimits::AddLimit( const EAxis pAxis,
53  const G4double pMin,
54  const G4double pMax )
55 {
56  if ( pAxis == kXAxis )
57  {
58  if ( pMin > fxAxisMin ) fxAxisMin = pMin ;
59  if ( pMax < fxAxisMax ) fxAxisMax = pMax ;
60  }
61  else if ( pAxis == kYAxis )
62  {
63  if ( pMin > fyAxisMin ) fyAxisMin = pMin ;
64  if ( pMax < fyAxisMax ) fyAxisMax = pMax ;
65  }
66  else
67  {
68  assert( pAxis == kZAxis ) ;
69 
70  if ( pMin > fzAxisMin ) fzAxisMin = pMin ;
71  if ( pMax < fzAxisMax ) fzAxisMax = pMax ;
72  }
73 }
74 
76 //
77 // ClipToLimits
78 //
79 // Clip the line segment pStart->pEnd to the volume described by the
80 // current limits. Return true if the line remains after clipping,
81 // else false, and leave the vectors in an undefined state.
82 //
83 // Process:
84 //
85 // Use Cohen-Sutherland clipping in 3D
86 // [Fundamentals of Interactive Computer Graphics,Foley & Van Dam]
87 //
89  G4ThreeVector& pEnd ) const
90 {
91  G4int sCode, eCode ;
92  G4bool remainsAfterClip ;
93 
94  // Determine if line is trivially inside (both outcodes==0) or outside
95  // (logical AND of outcodes !=0)
96 
97  sCode = OutCode(pStart) ;
98  eCode = OutCode(pEnd) ;
99 
100  if ( sCode & eCode )
101  {
102  // Trivially outside, no intersection with region
103 
104  remainsAfterClip = false;
105  }
106  else if ( sCode == 0 && eCode == 0 )
107  {
108  // Trivially inside, no intersections
109 
110  remainsAfterClip = true ;
111  }
112  else
113  {
114  // Line segment *may* cut volume boundaries
115  // At most, one end point is inside
116 
117  G4double x1, y1, z1, x2, y2, z2 ;
118 
119  x1 = pStart.x() ;
120  y1 = pStart.y() ;
121  z1 = pStart.z() ;
122 
123  x2 = pEnd.x() ;
124  y2 = pEnd.y() ;
125  z2 = pEnd.z() ;
126 
127  while ( sCode != eCode ) // Loop checking, 06.08.2015, G.Cosmo
128  {
129  // Copy vectors to work variables x1-z1,x2-z2
130  // Ensure x1-z1 lies outside volume, swapping vectors and outcodes
131  // if necessary
132 
133  if ( sCode )
134  {
135  if ( sCode & 0x01 ) // Clip against fxAxisMin
136  {
137  z1 += (fxAxisMin-x1)*(z2-z1)/(x2-x1);
138  y1 += (fxAxisMin-x1)*(y2-y1)/(x2-x1);
139  x1 = fxAxisMin;
140  }
141  else if ( sCode & 0x02 ) // Clip against fxAxisMax
142  {
143  z1 += (fxAxisMax-x1)*(z2-z1)/(x2-x1);
144  y1 += (fxAxisMax-x1)*(y2-y1)/(x2-x1);
145  x1 = fxAxisMax ;
146  }
147  else if ( sCode & 0x04 ) // Clip against fyAxisMin
148  {
149  x1 += (fyAxisMin-y1)*(x2-x1)/(y2-y1);
150  z1 += (fyAxisMin-y1)*(z2-z1)/(y2-y1);
151  y1 = fyAxisMin;
152  }
153  else if ( sCode & 0x08 ) // Clip against fyAxisMax
154  {
155  x1 += (fyAxisMax-y1)*(x2-x1)/(y2-y1);
156  z1 += (fyAxisMax-y1)*(z2-z1)/(y2-y1);
157  y1 = fyAxisMax;
158  }
159  else if ( sCode & 0x10 ) // Clip against fzAxisMin
160  {
161  x1 += (fzAxisMin-z1)*(x2-x1)/(z2-z1);
162  y1 += (fzAxisMin-z1)*(y2-y1)/(z2-z1);
163  z1 = fzAxisMin;
164  }
165  else if ( sCode & 0x20 ) // Clip against fzAxisMax
166  {
167  x1 += (fzAxisMax-z1)*(x2-x1)/(z2-z1);
168  y1 += (fzAxisMax-z1)*(y2-y1)/(z2-z1);
169  z1 = fzAxisMax;
170  }
171  }
172  if ( eCode ) // Clip 2nd end: repeat of 1st, but 1<>2
173  {
174  if ( eCode & 0x01 ) // Clip against fxAxisMin
175  {
176  z2 += (fxAxisMin-x2)*(z1-z2)/(x1-x2);
177  y2 += (fxAxisMin-x2)*(y1-y2)/(x1-x2);
178  x2 = fxAxisMin;
179  }
180  else if ( eCode & 0x02 ) // Clip against fxAxisMax
181  {
182  z2 += (fxAxisMax-x2)*(z1-z2)/(x1-x2);
183  y2 += (fxAxisMax-x2)*(y1-y2)/(x1-x2);
184  x2 = fxAxisMax;
185  }
186  else if ( eCode & 0x04 ) // Clip against fyAxisMin
187  {
188  x2 += (fyAxisMin-y2)*(x1-x2)/(y1-y2);
189  z2 += (fyAxisMin-y2)*(z1-z2)/(y1-y2);
190  y2 = fyAxisMin;
191  }
192  else if (eCode&0x08) // Clip against fyAxisMax
193  {
194  x2 += (fyAxisMax-y2)*(x1-x2)/(y1-y2);
195  z2 += (fyAxisMax-y2)*(z1-z2)/(y1-y2);
196  y2 = fyAxisMax;
197  }
198  else if ( eCode & 0x10 ) // Clip against fzAxisMin
199  {
200  x2 += (fzAxisMin-z2)*(x1-x2)/(z1-z2);
201  y2 += (fzAxisMin-z2)*(y1-y2)/(z1-z2);
202  z2 = fzAxisMin;
203  }
204  else if ( eCode & 0x20 ) // Clip against fzAxisMax
205  {
206  x2 += (fzAxisMax-z2)*(x1-x2)/(z1-z2);
207  y2 += (fzAxisMax-z2)*(y1-y2)/(z1-z2);
208  z2 = fzAxisMax;
209  }
210  }
211  pStart = G4ThreeVector(x1,y1,z1);
212  pEnd = G4ThreeVector(x2,y2,z2);
213  sCode = OutCode(pStart);
214  eCode = OutCode(pEnd);
215  }
216  if ( sCode == 0 && eCode == 0 ) remainsAfterClip = true;
217  else remainsAfterClip = false;
218  }
219  return remainsAfterClip;
220 }
221 
223 //
224 // Calculate the `outcode' for the specified vector:
225 // The following bits are set:
226 // 0 pVec.x()<fxAxisMin && IsXLimited()
227 // 1 pVec.x()>fxAxisMax && IsXLimited()
228 // 2 pVec.y()<fyAxisMin && IsYLimited()
229 // 3 pVec.y()>fyAxisMax && IsYLimited()
230 // 4 pVec.z()<fzAxisMin && IsZLimited()
231 // 5 pVec.z()>fzAxisMax && IsZLimited()
232 //
234 {
235  G4int code = 0 ; // The outcode
236 
237  if ( IsXLimited() )
238  {
239  if ( pVec.x() < fxAxisMin ) code |= 0x01 ;
240  if ( pVec.x() > fxAxisMax ) code |= 0x02 ;
241  }
242  if ( IsYLimited() )
243  {
244  if ( pVec.y() < fyAxisMin ) code |= 0x04 ;
245  if ( pVec.y() > fyAxisMax ) code |= 0x08 ;
246  }
247  if (IsZLimited())
248  {
249  if ( pVec.z() < fzAxisMin ) code |= 0x10 ;
250  if ( pVec.z() > fzAxisMax ) code |= 0x20 ;
251  }
252  return code;
253 }
254 
256 
257 std::ostream& operator << (std::ostream& os, const G4VoxelLimits& pLim)
258 {
259  os << "{";
260  if (pLim.IsXLimited())
261  {
262  os << "(" << pLim.GetMinXExtent()
263  << "," << pLim.GetMaxXExtent() << ") ";
264  }
265  else
266  {
267  os << "(-,-) ";
268  }
269  if (pLim.IsYLimited())
270  {
271  os << "(" << pLim.GetMinYExtent()
272  << "," << pLim.GetMaxYExtent() << ") ";
273  }
274  else
275  {
276  os << "(-,-) ";
277  }
278  if (pLim.IsZLimited())
279  {
280  os << "(" << pLim.GetMinZExtent()
281  << "," << pLim.GetMaxZExtent() << ")";
282  }
283  else
284  {
285  os << "(-,-)";
286  }
287  os << "}";
288  return os;
289 }