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 #include "G4BREPSolidBox.hh"
00037 #include "G4FPlane.hh"
00038 #include "G4Point3DVector.hh"
00039
00040 G4BREPSolidBox::G4BREPSolidBox(const G4String& name,
00041 const G4Point3D& Pt1,
00042 const G4Point3D& Pt2,
00043 const G4Point3D& Pt3,
00044 const G4Point3D& Pt4,
00045 const G4Point3D& Pt5,
00046 const G4Point3D& Pt6,
00047 const G4Point3D& Pt7,
00048 const G4Point3D& Pt8): G4BREPSolid(name)
00049 {
00050 nb_of_surfaces=6;
00051 active=1; PlaneSolid=1;
00052
00053
00054 constructorParams[0] = Pt1;
00055 constructorParams[1] = Pt2;
00056 constructorParams[2] = Pt3;
00057 constructorParams[3] = Pt4;
00058 constructorParams[4] = Pt5;
00059 constructorParams[5] = Pt6;
00060 constructorParams[6] = Pt7;
00061 constructorParams[7] = Pt8;
00062
00063 InitializeBox();
00064 }
00065
00066 G4BREPSolidBox::G4BREPSolidBox( __void__& a )
00067 : G4BREPSolid(a)
00068 {
00069 }
00070
00071 G4BREPSolidBox::~G4BREPSolidBox()
00072 {
00073 }
00074
00075 G4BREPSolidBox::G4BREPSolidBox(const G4BREPSolidBox& rhs)
00076 : G4BREPSolid(rhs), Rotation(rhs.Rotation)
00077 {
00078 for (size_t i=0; i<8; ++i) { constructorParams[i]= rhs.constructorParams[i]; }
00079 InitializeBox();
00080 }
00081
00082 G4BREPSolidBox& G4BREPSolidBox::operator = (const G4BREPSolidBox& rhs)
00083 {
00084
00085
00086 if (this == &rhs) { return *this; }
00087
00088
00089
00090 G4BREPSolid::operator=(rhs);
00091
00092
00093
00094 Rotation= rhs.Rotation;
00095 for (size_t i=0; i<8; ++i) { constructorParams[i]= rhs.constructorParams[i]; }
00096 InitializeBox();
00097
00098 return *this;
00099 }
00100
00101 void G4BREPSolidBox::InitializeBox()
00102 {
00103 SurfaceVec = new G4Surface*[6];
00104 G4Point3DVector PVec(4);
00105 G4int sense=0;
00106
00107 PVec[0] = constructorParams[0];
00108 PVec[1] = constructorParams[1];
00109 PVec[2] = constructorParams[2];
00110 PVec[3] = constructorParams[3];
00111 SurfaceVec[0] = new G4FPlane(&PVec);
00112
00113 PVec[2] = constructorParams[5];
00114 PVec[3] = constructorParams[4];
00115 SurfaceVec[1] = new G4FPlane(&PVec,0,sense);
00116
00117 PVec[0] = constructorParams[1];
00118 PVec[1] = constructorParams[5];
00119 PVec[2] = constructorParams[6];
00120 PVec[3] = constructorParams[2];
00121 SurfaceVec[2] = new G4FPlane(&PVec);
00122
00123 PVec[0] = constructorParams[2];
00124 PVec[1] = constructorParams[6];
00125 PVec[2] = constructorParams[7];
00126 PVec[3] = constructorParams[3];
00127 SurfaceVec[3] = new G4FPlane(&PVec);
00128
00129 PVec[0] = constructorParams[0];
00130 PVec[1] = constructorParams[4];
00131 PVec[2] = constructorParams[7];
00132 PVec[3] = constructorParams[3];
00133 SurfaceVec[4] = new G4FPlane(&PVec,0,sense);
00134
00135 PVec[0] = constructorParams[4];
00136 PVec[1] = constructorParams[5];
00137 PVec[2] = constructorParams[6];
00138 PVec[3] = constructorParams[7];
00139 SurfaceVec[5] = new G4FPlane(&PVec,0,sense);
00140
00141 Initialize();
00142 }
00143
00144 EInside G4BREPSolidBox::Inside(register const G4ThreeVector& Pt) const
00145 {
00146 G4Point3D Point(Pt);
00147
00148
00149 G4Point3D min = bbox->GetBoxMin();
00150 min += G4Point3D(-0.5*kCarTolerance,-0.5*kCarTolerance,-0.5*kCarTolerance);
00151
00152 G4Point3D max = bbox->GetBoxMax();
00153 max += G4Point3D(0.5*kCarTolerance,0.5*kCarTolerance,0.5*kCarTolerance);
00154
00155 if( (Point.x() < min.x() || Point.x() > max.x()) ||
00156 (Point.y() < min.y() || Point.y() > max.y()) ||
00157 (Point.z() < min.z() || Point.z() > max.z()) )
00158 return kOutside;
00159
00160 if( (Point.x() > min.x() && Point.x() < max.x())&&
00161 (Point.y() > min.y() && Point.y() < max.y())&&
00162 (Point.z() > min.z() && Point.z() < max.z()) )
00163 return kInside;
00164
00165 return kSurface;
00166 }
00167
00168 G4VSolid* G4BREPSolidBox::Clone() const
00169 {
00170 return new G4BREPSolidBox(*this);
00171 }
00172
00173 std::ostream& G4BREPSolidBox::StreamInfo(std::ostream& os) const
00174 {
00175
00176
00177 G4BREPSolid::StreamInfo( os )
00178 << "\n"
00179 << " Pt1: " << constructorParams[0]
00180 << " Pt2: " << constructorParams[1]
00181 << " Pt3: " << constructorParams[2]
00182 << " Pt4: " << constructorParams[3]
00183 << "\n Pt5: " << constructorParams[4]
00184 << " Pt6: " << constructorParams[5]
00185 << " Pt7: " << constructorParams[6]
00186 << " Pt8: " << constructorParams[7]
00187 << "\n-----------------------------------------------------------\n";
00188
00189 return os;
00190 }
00191