Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/examples/advanced/fastAerosol/include/FastAerosol.hh

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ] 

  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 // --------------------------------------------------------------------
 28 // GEANT 4 class header file
 29 //
 30 // 
 31 // FastAerosol
 32 //
 33 // Class description:
 34 //
 35 //   A FastAerosol is a collection of points in a voxelized 
 36 //   arbitrarily-shaped volume with methods implementing population of
 37 //   grids/voxels and for efficiently finding the nearest point.
 38 //
 39 // Author: A.Knaian (ara@nklabs.com), N.MacFadden (natemacfadden@gmail.com)
 40 // --------------------------------------------------------------------
 41 
 42 #ifndef FastAerosol_h
 43 #define FastAerosol_h
 44 
 45 #include "globals.hh"
 46 #include "Randomize.hh"
 47 #include "G4ThreeVector.hh"
 48 #include "G4VSolid.hh"
 49 
 50 // rotations and number density distribution
 51 #include <functional>
 52 #include "G4RotationMatrix.hh"
 53 #include <atomic>
 54 
 55 //using namespace std;
 56 
 57 class FastAerosol 
 58 {
 59 public:
 60  // Constructor; creates a random cloud of droplets
 61  FastAerosol(const G4String& pName, G4VSolid* pCloud,
 62        G4double pR, G4double pMinD, 
 63        G4double pAvgNumDens, G4double pdR,
 64        std::function<G4double (G4ThreeVector)> pNumDensDistribution);
 65 
 66 FastAerosol(const G4String& pName, G4VSolid* pCloud,
 67       G4double pR, G4double pMinD, 
 68       G4double pNumDens, G4double pdR);
 69 
 70 FastAerosol(const G4String& pName, G4VSolid* pCloud,
 71       G4double pR, G4double pMinD, G4double pNumDens);
 72     
 73 ~FastAerosol()=default;
 74 
 75 // Populate all grids. Otherwise, they are populated on-the-fly
 76 void PopulateAllGrids();
 77 
 78 // Save locations of droplets to a file for visualization/analysis purposes
 79 void SaveToFile(const char *filename);
 80 
 81 // Get absolutely nearest droplet - must be public as FastAerosolSolid uses it
 82 bool GetNearestDroplet(const G4ThreeVector &p, G4ThreeVector &center, G4double &closestDistance, G4double stepLim, G4VSolid* droplet, std::function<G4RotationMatrix (G4ThreeVector)> rotation);
 83     
 84 // Get nearest droplet along a vector - must be public as FastAerosolSolid uses it
 85 bool GetNearestDroplet(const G4ThreeVector &p, const G4ThreeVector &v, G4ThreeVector &center, G4double &closestDistance, G4double stepLim, G4VSolid* droplet, std::function<G4RotationMatrix (G4ThreeVector)> rotation);
 86 
 87 // ======
 88 // Inline
 89 // ======
 90 // Input quantities
 91 inline G4String GetName() const; //fasterosol name
 92 inline G4VSolid* GetBulk() const; // bulk shape
 93 inline G4double GetRadius() const; // droplet radius
 94 inline G4double GetAvgNumDens() const;  // droplet number density
 95 //inline G4double GetPitch() const;   // grid pitch
 96 
 97 inline G4int GetNumDroplets() const;
 98 // in case the absolute number is more relevant than density
 99 
100 // Bulk quantities
101 inline G4double GetXHalfLength() const;
102 inline G4double GetYHalfLength() const;
103 inline G4double GetZHalfLength() const;
104 inline void GetBoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const;
105 inline G4double GetCubicVolume() const;
106 
107 inline G4double DistanceToCloud(const G4ThreeVector &p);
108 inline G4double DistanceToCloud(const G4ThreeVector &p, const G4ThreeVector &v);
109 
110 // Misc getters and setters
111 inline long GetSeed();
112 inline void SetSeed(long seed);
113 
114 inline G4int GetNumPlacementTries();
115 inline void SetNumPlacementTries(G4int numTries);
116 
117 inline G4int GetPreSphereR();
118 inline void SetPreSphereR(G4int fPreSphereRIn);
119 
120 inline std::function<G4double (G4ThreeVector)> GetDistribution(); // the droplet number density distribution
121 
122 inline G4double GetDropletsPerVoxel();
123 inline void SetDropletsPerVoxel(G4double newDropletsPerVoxel);
124 
125 // Printing diagnostic tool
126 inline void PrintPopulationReport();
127     
128 private:
129  G4double kCarTolerance;
130 
131 // Parameters, set in constructor
132  G4String fName;
133  G4VSolid* fCloud; // Solid volume of the cloud
134  G4double fDx, fDy, fDz; // Half widths
135  G4double fR; // Bounding radius of each droplet
136  G4double fR2;  // Bounding radius squared of each droplet
137  G4double fdR;  // Uncertainty in DistanceToIn droplet when just using knowledge of droplet center
138 
139 G4double fMinD; // Minimum distance allowed between faces of droplets when constructing random array of droplets
140 
141 std::function<G4double (G4ThreeVector)> fDistribution;
142 G4double fAvgNumDens; // Average droplet number density
143 
144 long int fNumDroplets = 0;
145 // Number of droplets that have been created
146 
147 G4double fGridPitch;        
148 // Pitch of collision detection grid.  Must be greater than diameter of droplets for correctness of collision detection.
149 
150 // Ramdom engine
151 CLHEP::HepJamesRandom fCloudEngine;
152 long fSeed = 0; // Global random seed
153 
154 G4double fDropletsPerVoxel = 4.0; 
155 // Expected number of droplets per voxel
156 
157 // How far the voxel center must be inside the bulk 
158 //order for there to be no risk of placing a droplet outside
159 G4double fEdgeDistance;
160 
161 // Grid variables
162 std::vector<std::vector<G4ThreeVector>> fGrid;
163 // Grid of lists of inidices to grid points, 
164 //used for fast collsion checking
165 
166 std::vector<G4double> fGridMean;      
167 // Array listing mean count for each voxel
168 
169 std::atomic<bool> *fGridValid;    
170 // Array listing validity of each grid. uses atomic variables
171 
172 G4int fNx, fNy, fNz; // Number of x, y, and z elements in fGrid
173 
174 G4int fNxy; // Cached fNx*fNy
175 long int fNumGridCells; // Cached fNx*fNy*fNz
176 
177 G4double fCollisionLimit2;
178 // Threshold distance squared when checking for collsion
179     
180 G4int fNumNewPointTries = 100;    
181 // How many times we try to place droplets
182     
183 G4double fMaxDropPercent = 1.0;   
184 // The maximal percentage of skipped droplets before crashing0
185 
186 G4int fMaxDropCount;        
187 // The maximal number of skipped droplets before crashing
188     
189 G4int fNumDropped = 0;        
190 // Number of skipped droplets due to collisions/out of bulk placement
191 
192 G4int fNumCollisions = 0;     
193 // How many collisions occured when attempting to place
194 
195 // Droplet search variables
196 G4int fVectorSearchRadius;      
197 // maximum vector search radius
198 
199 // Droplet placement functions
200 // ===========================
201 void InitializeGrid();
202 
203 G4bool FindNewPoint(G4bool edgeVoxel, G4double dX, G4double dY, G4double dZ, G4double minX, G4double minY, G4double minZ, G4ThreeVector &foundVec);
204 
205 G4double VoxelOverlap(G4ThreeVector voxelCenter, G4int nStat, G4double epsilon);
206 
207 bool CheckCollision(G4double x, G4double y, G4double z);
208 bool CheckCollisionInsideGrid(G4double x, G4double y, G4double z, unsigned int xi, unsigned int yi, unsigned int zi);
209 bool CheckCollisionWithDroplet(G4double x, G4double y, G4double z, G4ThreeVector p);
210 
211 // Droplet distance functions
212 // ==========================
213 void SearchSphere(G4int searchRad, G4double &minDistance, std::vector<G4ThreeVector> &candidates, std::vector<G4double> &distances2, G4int xGrid, G4int yGrid, G4int zGrid, const G4ThreeVector &p);
214 
215 void GetNearestDropletInsideRegion(G4double &minDistance, G4ThreeVector &center, int xGrid, int yGrid, int zGrid, int xWidth, int yWidth, int zWidth, const G4ThreeVector &p, const G4ThreeVector &v, G4VSolid* droplet, std::function<G4RotationMatrix (G4ThreeVector)> rotation);
216 
217 void GetNearestDropletInsideGrid(G4double &minDistance, std::vector<G4ThreeVector> &candidates, std::vector<G4double> &distances2, unsigned int xGrid, unsigned int yGrid, unsigned int zGrid, const G4ThreeVector &p);
218 
219 void GetNearestDropletInsideGrid(G4double &minDistance, G4ThreeVector &center, unsigned int xGrid, unsigned int yGrid, unsigned int zGrid, const G4ThreeVector &p, const G4ThreeVector &v, G4VSolid* droplet, std::function<G4RotationMatrix (G4ThreeVector)> rotation);
220 
221 // Voxelized sphere methods
222 // ========================
223 // a collection of points as in {{x1,y1},{x2,y2},...}
224 typedef std::vector<std::vector<int>> fCircleType;
225 
226 // a collection of points describing a spherical shell
227 // with points (x,y,z)=(i-R,j-R,sphere[i][j][k])
228 // that is, first index gives x-position, second index
229 // gives y-position, and the value gives z-position
230 //
231 // this is done so that searching may be optimized:
232 // if searching some x=i-R that is outside the
233 // aerosol's bounding box, immediately increment i
234 // (similar for y).
235 typedef std::vector<std::vector<std::vector<int>>> fSphereType;
236 
237 G4int fMaxCircleR;
238 G4int fMaxSphereR;
239 G4int fPreSphereR = 20;
240 std::vector<fCircleType> fCircleCollection;
241 std::vector<fSphereType> fSphereCollection;
242 fSphereType MakeSphere(G4int R);
243 fCircleType MakeCircle(G4int R);
244 fCircleType MakeHalfCircle(G4int R);
245 
246 void PopulateGrid(unsigned int xi, unsigned int yi, unsigned int zi, unsigned int& gi); 
247 
248 // ======
249 // Inline
250 // ======
251 inline bool GetGrid(const G4ThreeVector &p, G4int &xGrid, G4int &yGrid, G4int &zGrid);
252 
253 inline bool AnyIndexOutOfBounds(G4int xGrid, G4int yGrid, G4int zGrid);
254 
255 inline unsigned int GetGridIndex(unsigned int xi, unsigned int yi, unsigned int zi);
256 
257 inline G4ThreeVector GetIndexCoord(G4int index);
258 
259 inline std::pair<G4int, G4int> GetMinMaxSide(G4int index, G4int numGrids);  
260 };
261 
262 #include "FastAerosol.icc"
263 
264 #endif
265