File indexing completed on 2025-01-18 09:58:55
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047
0048
0049
0050
0051
0052
0053
0054
0055
0056
0057
0058
0059
0060
0061 #ifndef G4PHYSICALVOLUMEMODEL_HH
0062 #define G4PHYSICALVOLUMEMODEL_HH
0063
0064 #include "G4VModel.hh"
0065 #include "G4ModelingParameters.hh"
0066
0067 #include "G4VTouchable.hh"
0068 #include "G4Transform3D.hh"
0069 #include "G4Plane3D.hh"
0070 #include <iostream>
0071 #include <vector>
0072 #include <map>
0073
0074 class G4VPhysicalVolume;
0075 class G4LogicalVolume;
0076 class G4VSolid;
0077 class G4Material;
0078 class G4VisAttributes;
0079 class G4AttDef;
0080 class G4AttValue;
0081
0082 class G4PhysicalVolumeModel: public G4VModel {
0083
0084 public:
0085
0086 enum {UNLIMITED = -1};
0087
0088 enum ClippingMode {subtraction, intersection};
0089
0090
0091 class G4PhysicalVolumeNodeID {
0092 public:
0093 G4PhysicalVolumeNodeID
0094 (G4VPhysicalVolume* pPV = 0,
0095 G4int iCopyNo = 0,
0096 G4int depth = 0,
0097 const G4Transform3D& transform = G4Transform3D(),
0098 G4bool drawn = true):
0099 fpPV(pPV),
0100 fCopyNo(iCopyNo),
0101 fNonCulledDepth(depth),
0102 fTransform(transform),
0103 fDrawn(drawn) {}
0104 G4VPhysicalVolume* GetPhysicalVolume() const {return fpPV;}
0105 G4int GetCopyNo() const {return fCopyNo;}
0106 G4int GetNonCulledDepth() const {return fNonCulledDepth;}
0107 const G4Transform3D& GetTransform() const {return fTransform;}
0108 G4bool GetDrawn() const {return fDrawn;}
0109 void SetPhysicalVolume(G4VPhysicalVolume* v) {fpPV = v;}
0110 void SetCopyNo (G4int n) {fCopyNo = n;}
0111 void SetNonCulledDepth(G4int d) {fNonCulledDepth = d;}
0112 void SetTransform (const G4Transform3D& t) {fTransform = t;}
0113 void SetDrawn (G4bool b) {fDrawn = b;}
0114 G4bool operator< (const G4PhysicalVolumeNodeID& right) const;
0115 G4bool operator!= (const G4PhysicalVolumeNodeID& right) const;
0116 G4bool operator== (const G4PhysicalVolumeNodeID& right) const {
0117 return !operator!= (right);
0118 }
0119 private:
0120 G4VPhysicalVolume* fpPV;
0121 G4int fCopyNo;
0122 G4int fNonCulledDepth;
0123 G4Transform3D fTransform;
0124 G4bool fDrawn;
0125 };
0126
0127
0128 class G4PhysicalVolumeModelTouchable: public G4VTouchable {
0129 public:
0130 G4PhysicalVolumeModelTouchable
0131 (const std::vector<G4PhysicalVolumeNodeID>& fullPVPath);
0132 const G4ThreeVector& GetTranslation(G4int depth) const;
0133 const G4RotationMatrix* GetRotation(G4int depth) const;
0134 G4VPhysicalVolume* GetVolume(G4int depth) const;
0135 G4VSolid* GetSolid(G4int depth) const;
0136 G4int GetReplicaNumber(G4int depth) const;
0137 G4int GetHistoryDepth() const {return G4int(fFullPVPath.size());}
0138 private:
0139 const std::vector<G4PhysicalVolumeNodeID>& fFullPVPath;
0140 };
0141
0142
0143 struct TouchableProperties {
0144 TouchableProperties(): fpTouchablePV(nullptr), fCopyNo(0) {}
0145 G4ModelingParameters::PVNameCopyNoPath fTouchablePath;
0146 G4VPhysicalVolume* fpTouchablePV;
0147 G4int fCopyNo;
0148 G4Transform3D fTouchableGlobalTransform;
0149 std::vector<G4PhysicalVolumeNodeID> fTouchableBaseFullPVPath;
0150 std::vector<G4PhysicalVolumeNodeID> fTouchableFullPVPath;
0151 };
0152
0153 G4PhysicalVolumeModel
0154 (G4VPhysicalVolume* = 0,
0155 G4int requestedDepth = UNLIMITED,
0156 const G4Transform3D& modelTransformation = G4Transform3D(),
0157 const G4ModelingParameters* = 0,
0158 G4bool useFullExtent = false,
0159 const std::vector<G4PhysicalVolumeNodeID>& baseFullPVPath =
0160 std::vector<G4PhysicalVolumeNodeID>());
0161
0162 virtual ~G4PhysicalVolumeModel ();
0163
0164 void DescribeYourselfTo (G4VGraphicsScene&);
0165
0166
0167
0168
0169 G4String GetCurrentDescription () const;
0170
0171
0172 G4String GetCurrentTag () const;
0173
0174
0175 G4VPhysicalVolume* GetTopPhysicalVolume () const {return fpTopPV;}
0176
0177 G4int GetRequestedDepth () const {return fRequestedDepth;}
0178
0179 const G4VSolid* GetClippingSolid () const
0180 {return fpClippingSolid;}
0181
0182 G4int GetCurrentDepth() const {return fCurrentDepth;}
0183
0184
0185 const G4Transform3D& GetTransformation() const {return fTransform;}
0186
0187
0188 G4VPhysicalVolume* GetCurrentPV() const {return fpCurrentPV;}
0189
0190
0191 G4int GetCurrentPVCopyNo() const {return fCurrentPVCopyNo;}
0192
0193
0194 G4LogicalVolume* GetCurrentLV() const {return fpCurrentLV;}
0195
0196
0197 G4Material* GetCurrentMaterial() const {return fpCurrentMaterial;}
0198
0199
0200 const G4Transform3D& GetCurrentTransform() const {return fCurrentTransform;}
0201
0202
0203 const std::vector<G4PhysicalVolumeNodeID>& GetBaseFullPVPath() const
0204 {return fBaseFullPVPath;}
0205
0206
0207
0208
0209
0210 const std::vector<G4PhysicalVolumeNodeID>& GetFullPVPath() const
0211 {return fFullPVPath;}
0212
0213
0214
0215
0216
0217 const std::vector<G4PhysicalVolumeNodeID>& GetDrawnPVPath() const
0218 {return fDrawnPVPath;}
0219
0220
0221
0222
0223
0224 static G4ModelingParameters::PVNameCopyNoPath GetPVNameCopyNoPath
0225 (const std::vector<G4PhysicalVolumeNodeID>&);
0226
0227
0228 static G4String GetPVNamePathString(const std::vector<G4PhysicalVolumeNodeID>&);
0229
0230
0231
0232 const std::map<G4String,G4AttDef>* GetAttDefs() const;
0233
0234
0235 std::vector<G4AttValue>* CreateCurrentAttValues() const;
0236
0237
0238
0239
0240
0241
0242
0243
0244
0245 const std::map<G4int,G4int>& GetNumberOfTouchables() const {return fNTouchables;}
0246
0247
0248 void SetRequestedDepth (G4int requestedDepth) {
0249 fRequestedDepth = requestedDepth;
0250 }
0251
0252 void SetClippingSolid (G4VSolid* pClippingSolid) {
0253 fpClippingSolid = pClippingSolid;
0254 }
0255
0256 void SetClippingMode (ClippingMode mode) {
0257 fClippingMode = mode;
0258 }
0259
0260 G4bool Validate (G4bool warn);
0261
0262
0263 void Abort () const {fAbort = true;}
0264
0265
0266 void CurtailDescent() const {fCurtailDescent = true;}
0267
0268
0269 void CalculateExtent ();
0270
0271 protected:
0272
0273 void VisitGeometryAndGetVisReps (G4VPhysicalVolume*,
0274 G4int requestedDepth,
0275 const G4Transform3D&,
0276 G4VGraphicsScene&);
0277
0278 void DescribeAndDescend (G4VPhysicalVolume*,
0279 G4int requestedDepth,
0280 G4LogicalVolume*,
0281 G4VSolid*,
0282 G4Material*,
0283 const G4Transform3D&,
0284 G4VGraphicsScene&);
0285
0286 virtual void DescribeSolid (const G4Transform3D& theAT,
0287 G4VSolid* pSol,
0288 const G4VisAttributes* pVisAttribs,
0289 G4VGraphicsScene& sceneHandler);
0290
0291
0292
0293
0294 G4VPhysicalVolume* fpTopPV;
0295 G4String fTopPVName;
0296 G4int fTopPVCopyNo;
0297 G4int fRequestedDepth;
0298
0299 G4bool fUseFullExtent;
0300 G4Transform3D fTransform;
0301 G4int fCurrentDepth;
0302 G4VPhysicalVolume* fpCurrentPV;
0303 G4int fCurrentPVCopyNo;
0304 G4LogicalVolume* fpCurrentLV;
0305 G4Material* fpCurrentMaterial;
0306 G4Transform3D fCurrentTransform;
0307 std::vector<G4PhysicalVolumeNodeID> fBaseFullPVPath;
0308 std::vector<G4PhysicalVolumeNodeID> fFullPVPath;
0309 std::vector<G4PhysicalVolumeNodeID> fDrawnPVPath;
0310 mutable G4bool fAbort;
0311 mutable G4bool fCurtailDescent;
0312 G4VSolid* fpClippingSolid;
0313 ClippingMode fClippingMode;
0314 G4int fNClippers;
0315 std::map<G4int,G4int> fNTouchables;
0316
0317 private:
0318
0319
0320
0321 G4PhysicalVolumeModel (const G4PhysicalVolumeModel&);
0322 G4PhysicalVolumeModel& operator = (const G4PhysicalVolumeModel&);
0323 };
0324
0325 std::ostream& operator<<
0326 (std::ostream& os, const G4PhysicalVolumeModel::G4PhysicalVolumeNodeID&);
0327
0328 std::ostream& operator<<
0329 (std::ostream& os, const std::vector<G4PhysicalVolumeModel::G4PhysicalVolumeNodeID>&);
0330
0331 #endif