00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041 #ifndef G4SurfaceVoxelizer_HH
00042 #define G4SurfaceVoxelizer_HH
00043
00044 #include <vector>
00045 #include <string>
00046 #include <map>
00047
00048 #include "G4SurfBits.hh"
00049 #include "G4Box.hh"
00050 #include "G4VFacet.hh"
00051
00052 #include "Randomize.hh"
00053
00054 struct G4VoxelBox
00055 {
00056 G4ThreeVector hlen;
00057 G4ThreeVector pos;
00058 };
00059
00060 struct G4VoxelInfo
00061 {
00062 G4int count;
00063 G4int previous;
00064 G4int next;
00065 };
00066
00067 class G4SurfaceVoxelizer
00068 {
00069 friend class G4VoxelCandidatesIterator;
00070
00071 public:
00072
00073 template <typename T>
00074 static inline G4int BinarySearch(const std::vector<T> &vec, T value);
00075
00076 void Voxelize(std::vector<G4VFacet *> &facets);
00077
00078 void DisplayVoxelLimits();
00079 void DisplayBoundaries();
00080 void DisplayListNodes();
00081
00082 G4SurfaceVoxelizer();
00083 ~G4SurfaceVoxelizer();
00084
00085 void GetCandidatesVoxel(std::vector<G4int> &voxels);
00086
00087
00088
00089 G4int GetCandidatesVoxelArray(const G4ThreeVector &point,
00090 std::vector<G4int> &list,
00091 G4SurfBits *crossed=0) const;
00092
00093
00094 G4int GetCandidatesVoxelArray(const std::vector<G4int> &voxels,
00095 const G4SurfBits bitmasks[],
00096 std::vector<G4int> &list,
00097 G4SurfBits *crossed=0) const;
00098 G4int GetCandidatesVoxelArray(const std::vector<G4int> &voxels,
00099 std::vector<G4int> &list,
00100 G4SurfBits *crossed=0)const;
00101
00102 inline const std::vector<G4VoxelBox> &GetBoxes() const;
00103
00104
00105
00106 inline const std::vector<G4double> &GetBoundary(G4int index) const;
00107
00108 G4bool UpdateCurrentVoxel(const G4ThreeVector &point,
00109 const G4ThreeVector &direction,
00110 std::vector<G4int> &curVoxel) const;
00111
00112 inline void GetVoxel(std::vector<G4int> &curVoxel,
00113 const G4ThreeVector &point) const;
00114
00115 inline G4int GetBitsPerSlice () const;
00116
00117 G4bool Contains(const G4ThreeVector &point) const;
00118
00119 G4double DistanceToNext(const G4ThreeVector &point,
00120 const G4ThreeVector &direction,
00121 const std::vector<G4int> &curVoxel) const;
00122
00123 G4double DistanceToFirst(const G4ThreeVector &point,
00124 const G4ThreeVector &direction) const;
00125
00126 G4double DistanceToBoundingBox(const G4ThreeVector &point) const;
00127
00128 inline G4int GetVoxelsIndex(G4int x, G4int y, G4int z) const;
00129 inline G4int GetVoxelsIndex(const std::vector<G4int> &voxels) const;
00130
00131 inline G4int GetPointIndex(const G4ThreeVector &p) const;
00132
00133 inline const G4SurfBits &Empty() const;
00134 inline G4bool IsEmpty(G4int index) const;
00135
00136 void SetMaxVoxels(G4int max);
00137 void SetMaxVoxels(const G4ThreeVector &reductionRatio);
00138
00139 inline G4int GetMaxVoxels(G4ThreeVector &ratioOfReduction);
00140
00141 G4int AllocatedMemory();
00142
00143 inline long long GetCountOfVoxels() const;
00144
00145 inline long long CountVoxels(std::vector<G4double> boundaries[]) const;
00146
00147 inline G4int GetCandidates(std::vector<G4int> &curVoxel,
00148 std::vector<G4int> *&candidates,
00149 std::vector<G4int> &space) const;
00150 inline const std::vector<G4int> &
00151 GetCandidates(std::vector<G4int> &curVoxel) const;
00152
00153 inline G4int GetVoxelBoxesSize() const;
00154
00155 inline const G4VoxelBox &GetVoxelBox(G4int i) const;
00156
00157 inline const std::vector<G4int> &GetVoxelBoxCandidates(G4int i) const;
00158
00159 static G4double MinDistanceToBox (const G4ThreeVector &aPoint,
00160 const G4ThreeVector &f);
00161
00162 static G4int SetDefaultVoxelsCount(G4int count);
00163
00164 static G4int GetDefaultVoxelsCount();
00165
00166 private:
00167
00168 class G4VoxelComparator
00169 {
00170 public:
00171
00172 std::vector<G4VoxelInfo> &fVoxels;
00173
00174 G4VoxelComparator(std::vector<G4VoxelInfo> &voxels) : fVoxels(voxels) {}
00175
00176 G4bool operator()(const G4int& l, const G4int& r) const
00177 {
00178 G4VoxelInfo &lv = fVoxels[l], &rv = fVoxels[r];
00179 G4int left = lv.count + fVoxels[lv.next].count;
00180 G4int right = rv.count + fVoxels[rv.next].count;
00181 return (left == right) ? l < r : left < right;
00182 }
00183 };
00184
00185 private:
00186
00187 static int fDefaultVoxelsCount;
00188
00189 void BuildEmpty ();
00190
00191 G4String GetCandidatesAsString(const G4SurfBits &bits);
00192
00193 void CreateSortedBoundary(std::vector<G4double> &boundaryRaw, G4int axis);
00194
00195 void BuildBoundaries();
00196
00197 void BuildReduceVoxels(std::vector<G4double> fBoundaries[],
00198 G4ThreeVector reductionRatio);
00199 void BuildReduceVoxels2(std::vector<G4double> fBoundaries[],
00200 G4ThreeVector reductionRatio);
00201
00202 void BuildVoxelLimits(std::vector<G4VFacet *> &facets);
00203
00204 void DisplayBoundaries(std::vector<G4double> &fBoundaries);
00205
00206 void BuildBitmasks(std::vector<G4double> fBoundaries[],
00207 G4SurfBits bitmasks[], G4bool countsOnly=false);
00208
00209 void BuildBoundingBox();
00210
00211 void SetReductionRatio(G4int maxVoxels, G4ThreeVector &reductionRatio);
00212
00213 void CreateMiniVoxels(std::vector<G4double> fBoundaries[],
00214 G4SurfBits bitmasks[]);
00215 static void FindComponentsFastest(unsigned int mask,
00216 std::vector<G4int> &list, G4int i);
00217
00218 private:
00219
00220 std::vector<G4VoxelBox> fVoxelBoxes;
00221 std::vector<std::vector<G4int> > fVoxelBoxesCandidates;
00222 mutable std::map<G4int, std::vector<G4int> > fCandidates;
00223
00224 const std::vector<G4int> fNoCandidates;
00225
00226 long long fCountOfVoxels;
00227
00228 G4int fNPerSlice;
00229
00230 std::vector<G4VoxelBox> fBoxes;
00231
00232
00233 std::vector<G4double> fBoundaries[3];
00234
00235
00236 std::vector<G4int> fCandidatesCounts[3];
00237
00238 G4int fTotalCandidates;
00239
00240 G4SurfBits fBitmasks[3];
00241
00242 G4ThreeVector fBoundingBoxCenter;
00243
00244 G4Box fBoundingBox;
00245
00246 G4ThreeVector fBoundingBoxSize;
00247
00248 G4ThreeVector fReductionRatio;
00249
00250 G4int fMaxVoxels;
00251
00252 G4double fTolerance;
00253
00254 G4SurfBits fEmpty;
00255 };
00256
00257 #include "G4SurfaceVoxelizer.icc"
00258
00259 #endif