1 // Copyright (c) 1995-1999 Matra Datavision
2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
3 //
4 // This file is part of Open CASCADE Technology software library.
5 //
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License version 2.1 as published
8 // by the Free Software Foundation, with special exception defined in the file
9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10 // distribution for complete text of the license and disclaimer of any warranty.
11 //
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
14
15
16 #include <BRepAdaptor_Curve.hxx>
17 #include <BRepGProp_Cinert.hxx>
18 #include <BRepGProp_EdgeTool.hxx>
19 #include <gp_Pnt.hxx>
20 #include <math.hxx>
21 #include <TColStd_Array1OfReal.hxx>
22
BRepGProp_Cinert()23 BRepGProp_Cinert::BRepGProp_Cinert(){}
24
SetLocation(const gp_Pnt & CLocation)25 void BRepGProp_Cinert::SetLocation(const gp_Pnt& CLocation)
26 {
27 loc = CLocation;
28 }
29
Perform(const BRepAdaptor_Curve & C)30 void BRepGProp_Cinert::Perform (const BRepAdaptor_Curve& C)
31 {
32
33 Standard_Real Ix, Iy, Iz, Ixx, Iyy, Izz, Ixy, Ixz, Iyz;
34 dim = Ix = Iy = Iz = Ixx = Iyy = Izz = Ixy = Ixz = Iyz = 0.0;
35
36 Standard_Real Lower = BRepGProp_EdgeTool::FirstParameter (C);
37 Standard_Real Upper = BRepGProp_EdgeTool::LastParameter (C);
38 Standard_Integer Order = Min(BRepGProp_EdgeTool::IntegrationOrder (C),
39 math::GaussPointsMax());
40
41 gp_Pnt P; //value on the curve
42 gp_Vec V1; //first derivative on the curve
43 Standard_Real ds; //curvilign abscissae
44 Standard_Real ur, um, u;
45 Standard_Real x, y, z;
46 Standard_Real xloc, yloc, zloc;
47
48 math_Vector GaussP (1, Order);
49 math_Vector GaussW (1, Order);
50
51 //Recuperation des points de Gauss dans le fichier GaussPoints.
52 math::GaussPoints (Order,GaussP);
53 math::GaussWeights (Order,GaussW);
54
55 // modified by NIZHNY-MKK Thu Jun 9 12:13:21 2005.BEGIN
56 Standard_Integer nbIntervals = BRepGProp_EdgeTool::NbIntervals(C, GeomAbs_CN);
57 Standard_Boolean bHasIntervals = (nbIntervals > 1);
58 TColStd_Array1OfReal TI(1, nbIntervals + 1);
59
60 if(bHasIntervals) {
61 BRepGProp_EdgeTool::Intervals(C, TI, GeomAbs_CN);
62 }
63 else {
64 nbIntervals = 1;
65 }
66 Standard_Integer nIndex = 0;
67 Standard_Real UU1 = Min(Lower, Upper);
68 Standard_Real UU2 = Max(Lower, Upper);
69
70 for(nIndex = 1; nIndex <= nbIntervals; nIndex++) {
71 if(bHasIntervals) {
72 Lower = Max(TI(nIndex), UU1);
73 Upper = Min(TI(nIndex+1), UU2);
74 }
75 else {
76 Lower = UU1;
77 Upper = UU2;
78 }
79
80 Standard_Real dimLocal, IxLocal, IyLocal, IzLocal, IxxLocal, IyyLocal, IzzLocal, IxyLocal, IxzLocal, IyzLocal;
81 dimLocal = IxLocal = IyLocal = IzLocal = IxxLocal = IyyLocal = IzzLocal = IxyLocal = IxzLocal = IyzLocal = 0.0;
82 // modified by NIZHNY-MKK Thu Jun 9 12:13:32 2005.END
83
84 loc.Coord (xloc, yloc, zloc);
85
86 Standard_Integer i;
87
88 // Calcul des integrales aux points de gauss :
89 um = 0.5 * (Upper + Lower);
90 ur = 0.5 * (Upper - Lower);
91
92 for (i = 1; i <= Order; i++) {
93 u = um + ur * GaussP (i);
94 BRepGProp_EdgeTool::D1 (C,u, P, V1);
95 ds = V1.Magnitude();
96 P.Coord (x, y, z);
97 x -= xloc;
98 y -= yloc;
99 z -= zloc;
100 ds *= GaussW (i);
101 dimLocal += ds;
102 IxLocal += x * ds;
103 IyLocal += y * ds;
104 IzLocal += z * ds;
105 IxyLocal += x * y * ds;
106 IyzLocal += y * z * ds;
107 IxzLocal += x * z * ds;
108 x *= x;
109 y *= y;
110 z *= z;
111 IxxLocal += (y + z) * ds;
112 IyyLocal += (x + z) * ds;
113 IzzLocal += (x + y) * ds;
114 }
115 // modified by NIZHNY-MKK Thu Jun 9 12:13:47 2005.BEGIN
116 dimLocal *= ur;
117 IxLocal *= ur;
118 IyLocal *= ur;
119 IzLocal *= ur;
120 IxxLocal *= ur;
121 IyyLocal *= ur;
122 IzzLocal *= ur;
123 IxyLocal *= ur;
124 IxzLocal *= ur;
125 IyzLocal *= ur;
126
127 dim += dimLocal;
128 Ix += IxLocal;
129 Iy += IyLocal;
130 Iz += IzLocal;
131 Ixx += IxxLocal;
132 Iyy += IyyLocal;
133 Izz += IzzLocal;
134 Ixy += IxyLocal;
135 Ixz += IxzLocal;
136 Iyz += IyzLocal;
137 }
138 // modified by NIZHNY-MKK Thu Jun 9 12:13:55 2005.END
139
140 inertia = gp_Mat (gp_XYZ (Ixx, -Ixy, -Ixz),
141 gp_XYZ (-Ixy, Iyy, -Iyz),
142 gp_XYZ (-Ixz, -Iyz, Izz));
143
144 if (Abs(dim) < gp::Resolution())
145 g = P;
146 else
147 g.SetCoord (Ix/dim, Iy/dim, Iz/dim);
148 }
149
150
BRepGProp_Cinert(const BRepAdaptor_Curve & C,const gp_Pnt & CLocation)151 BRepGProp_Cinert::BRepGProp_Cinert (const BRepAdaptor_Curve& C,
152 const gp_Pnt& CLocation)
153 {
154 SetLocation(CLocation);
155 Perform(C);
156 }
157