0023892: Missing intersection edge between two faces.
[occt.git] / src / GProp / GProp_VGProps.gxx
CommitLineData
b311480e 1// Copyright (c) 1995-1999 Matra Datavision
973c2be1 2// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 3//
973c2be1 4// This file is part of Open CASCADE Technology software library.
b311480e 5//
d5f74e42 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
973c2be1 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.
b311480e 11//
973c2be1 12// Alternatively, this file may be used under the terms of Open CASCADE
13// commercial license or contractual agreement.
b311480e 14
7fd59977 15#include <Standard_NotImplemented.hxx>
16#include <math_Vector.hxx>
17#include <math.hxx>
18#include <gp_Pnt2d.hxx>
19#include <gp_Vec2d.hxx>
20#include <gp_Pnt.hxx>
21#include <gp_Vec.hxx>
22
23#include <TColStd_Array1OfReal.hxx>
24#include <Precision.hxx>
25class HMath_Vector{
26 math_Vector *pvec;
27 void operator=(const math_Vector&){}
28 public:
29 HMath_Vector(){ pvec = 0;}
30 HMath_Vector(math_Vector* pv){ pvec = pv;}
31 ~HMath_Vector(){ if(pvec != 0) delete pvec;}
32 void operator=(math_Vector* pv){ if(pvec != pv && pvec != 0) delete pvec; pvec = pv;}
33 Standard_Real& operator()(Standard_Integer i){ return (*pvec).operator()(i);}
34 const Standard_Real& operator()(Standard_Integer i) const{ return (*pvec).operator()(i);}
35 const math_Vector* operator->() const{ return pvec;}
36 math_Vector* operator->(){ return pvec;}
37 math_Vector* Init(Standard_Real v, Standard_Integer i = 0, Standard_Integer iEnd = 0){
38 if(pvec == 0) return pvec;
39 if(iEnd - i == 0) pvec->Init(v);
40 else for(; i <= iEnd; i++) pvec->operator()(i) = v;
41 return pvec;
42 }
43};
44
45//Minimal value of interval's range for computation | minimal value of "dim" | ...
46static Standard_Real EPS_PARAM = Precision::Angular(), EPS_DIM = 1.E-30, ERROR_ALGEBR_RATIO = 2.0/3.0;
47//Maximum of GaussPoints on a subinterval and maximum of subintervals
48static Standard_Integer GPM = math::GaussPointsMax(), SUBS_POWER = 32, SM = SUBS_POWER*GPM + 1;
49static Standard_Boolean IS_MIN_DIM = 1; // if the value equal 0 error of algorithm calculted by static moments
50
51static math_Vector LGaussP0(1,GPM), LGaussW0(1,GPM),
52 LGaussP1(1,RealToInt(Ceiling(ERROR_ALGEBR_RATIO*GPM))), LGaussW1(1,RealToInt(Ceiling(ERROR_ALGEBR_RATIO*GPM)));
53static HMath_Vector L1 = new math_Vector(1,SM), L2 = new math_Vector(1,SM),
54 DimL = new math_Vector(1,SM), ErrL = new math_Vector(1,SM), ErrUL = new math_Vector(1,SM,0.0),
55 IxL = new math_Vector(1,SM), IyL = new math_Vector(1,SM), IzL = new math_Vector(1,SM),
56 IxxL = new math_Vector(1,SM), IyyL = new math_Vector(1,SM), IzzL = new math_Vector(1,SM),
57 IxyL = new math_Vector(1,SM), IxzL = new math_Vector(1,SM), IyzL = new math_Vector(1,SM);
58
59static math_Vector* LGaussP[] = {&LGaussP0,&LGaussP1};
60static math_Vector* LGaussW[] = {&LGaussW0,&LGaussW1};
61
62static math_Vector UGaussP0(1,GPM), UGaussW0(1,GPM),
63 UGaussP1(1,RealToInt(Ceiling(ERROR_ALGEBR_RATIO*GPM))), UGaussW1(1,RealToInt(Ceiling(ERROR_ALGEBR_RATIO*GPM)));
64static HMath_Vector U1 = new math_Vector(1,SM), U2 = new math_Vector(1,SM),
65 DimU = new math_Vector(1,SM), ErrU = new math_Vector(1,SM,0.0),
66 IxU = new math_Vector(1,SM), IyU = new math_Vector(1,SM), IzU = new math_Vector(1,SM),
67 IxxU = new math_Vector(1,SM), IyyU = new math_Vector(1,SM), IzzU = new math_Vector(1,SM),
68 IxyU = new math_Vector(1,SM), IxzU = new math_Vector(1,SM), IyzU = new math_Vector(1,SM);
69
70static math_Vector* UGaussP[] = {&UGaussP0,&UGaussP1};
71static math_Vector* UGaussW[] = {&UGaussW0,&UGaussW1};
72
73static Standard_Integer FillIntervalBounds(Standard_Real A, Standard_Real B, const TColStd_Array1OfReal& Knots,
74 HMath_Vector& VA, HMath_Vector& VB)
75{
76 Standard_Integer i = 1, iEnd = Knots.Upper(), j = 1, k = 1;
77 VA(j++) = A;
78 for(; i <= iEnd; i++){
79 Standard_Real kn = Knots(i);
80 if(A < kn)
258ff83b 81 {
82 if(kn < B)
83 {
84 VA(j++) = VB(k++) = kn;
85 }
86 else
87 {
88 break;
89 }
90 }
7fd59977 91 }
92 VB(k) = B;
93 return k;
94}
95
96static inline Standard_Integer MaxSubs(Standard_Integer n, Standard_Integer coeff = SUBS_POWER){
97 return n = IntegerLast()/coeff < n? IntegerLast(): n*coeff + 1;
98}
99
100static Standard_Integer LFillIntervalBounds(Standard_Real A, Standard_Real B, const TColStd_Array1OfReal& Knots,
101 const Standard_Integer NumSubs)
102{
103 Standard_Integer iEnd = Knots.Upper(), jEnd = L1->Upper();
104
105// Modified by Sergey KHROMOV - Wed Mar 26 11:22:50 2003
106 iEnd = Max(iEnd, MaxSubs(iEnd-1,NumSubs));
107 if(iEnd - 1 > jEnd){
108// iEnd = MaxSubs(iEnd-1,NumSubs);
109// Modified by Sergey KHROMOV - Wed Mar 26 11:22:51 2003
110 L1 = new math_Vector(1,iEnd); L2 = new math_Vector(1,iEnd);
111 DimL = new math_Vector(1,iEnd); ErrL = new math_Vector(1,iEnd,0.0); ErrUL = new math_Vector(1,iEnd,0.0);
112 IxL = new math_Vector(1,iEnd); IyL = new math_Vector(1,iEnd); IzL = new math_Vector(1,iEnd);
113 IxxL = new math_Vector(1,iEnd); IyyL = new math_Vector(1,iEnd); IzzL = new math_Vector(1,iEnd);
114 IxyL = new math_Vector(1,iEnd); IxzL = new math_Vector(1,iEnd); IyzL = new math_Vector(1,iEnd);
115 }
116 return FillIntervalBounds(A, B, Knots, L1, L2);
117}
118
119static Standard_Integer UFillIntervalBounds(Standard_Real A, Standard_Real B, const TColStd_Array1OfReal& Knots,
120 const Standard_Integer NumSubs)
121{
122 Standard_Integer iEnd = Knots.Upper(), jEnd = U1->Upper();
123
124// Modified by Sergey KHROMOV - Wed Mar 26 11:22:50 2003
125 iEnd = Max(iEnd, MaxSubs(iEnd-1,NumSubs));
126 if(iEnd - 1 > jEnd){
127// iEnd = MaxSubs(iEnd-1,NumSubs);
128// Modified by Sergey KHROMOV - Wed Mar 26 11:22:51 2003
129 U1 = new math_Vector(1,iEnd); U2 = new math_Vector(1,iEnd);
130 DimU = new math_Vector(1,iEnd); ErrU = new math_Vector(1,iEnd,0.0);
131 IxU = new math_Vector(1,iEnd); IyU = new math_Vector(1,iEnd); IzU = new math_Vector(1,iEnd);
132 IxxU = new math_Vector(1,iEnd); IyyU = new math_Vector(1,iEnd); IzzU = new math_Vector(1,iEnd);
133 IxyU = new math_Vector(1,iEnd); IxzU = new math_Vector(1,iEnd); IyzU = new math_Vector(1,iEnd);
134 }
135 return FillIntervalBounds(A, B, Knots, U1, U2);
136}
137
138static Standard_Real CCompute(Face& S, Domain& D, const Standard_Boolean ByPoint, const Standard_Real Coeff[],
139 const gp_Pnt& loc, Standard_Real& Dim, gp_Pnt& g, gp_Mat& inertia,
140 const Standard_Real EpsDim,
141 const Standard_Boolean isErrorCalculation, const Standard_Boolean isVerifyComputation)
142{
143 Standard_Boolean isNaturalRestriction = S.NaturalRestriction();
144
145 Standard_Integer NumSubs = SUBS_POWER;
146 Standard_Boolean isMinDim = IS_MIN_DIM;
147
148 Standard_Real Ix, Iy, Iz, Ixx, Iyy, Izz, Ixy, Ixz, Iyz;
149 Dim = Ix = Iy = Iz = Ixx = Iyy = Izz = Ixy = Ixz = Iyz = 0.0;
150 //boundary curve parametrization
151 Standard_Real l1, l2, lm, lr, l;
152 //Face parametrization in U and V direction
153 Standard_Real BV1, BV2, v;
154 Standard_Real BU1, BU2, u1, u2, um, ur, u;
155 S.Bounds (BU1, BU2, BV1, BV2); u1 = BU1;
156 //location point used to compute the inertia
157 Standard_Real xloc, yloc, zloc;
158 loc.Coord (xloc, yloc, zloc);
159 //location point used to compute the inertiard (xloc, yloc, zloc);
160 //Jacobien (x, y, z) -> (u, v) = ||n||
161 Standard_Real xn, yn, zn, s, ds, dDim;
162 Standard_Real x, y, z, xi, px, py, pz, yi, zi, d1, d2, d3;
163 //On the Face
164 gp_Pnt Ps;
165 gp_Vec VNor;
166 //On the boundary curve u-v
167 gp_Pnt2d Puv;
168 gp_Vec2d Vuv;
169 Standard_Real Dul; // Dul = Du / Dl
170 Standard_Real CDim[2], CIx, CIy, CIz, CIxx[2], CIyy[2], CIzz[2], CIxy, CIxz, CIyz;
171 Standard_Real LocDim[2], LocIx[2], LocIy[2], LocIz[2], LocIxx[2], LocIyy[2], LocIzz[2], LocIxy[2], LocIxz[2], LocIyz[2];
172
173 Standard_Integer iD = 0, NbLSubs, iLS, iLSubEnd, iGL, iGLEnd, NbLGaussP[2], LRange[2], iL, kL, kLEnd, IL, JL;
174 Standard_Integer i, NbUSubs, iUS, iUSubEnd, iGU, iGUEnd, NbUGaussP[2], URange[2], iU, kU, kUEnd, IU, JU;
175 Standard_Integer UMaxSubs, LMaxSubs;
176
177 Standard_Real ErrorU, ErrorL, ErrorLMax = 0.0, Eps=0.0, EpsL=0.0, EpsU=0.0;
178 iGLEnd = isErrorCalculation? 2: 1;
179
180 for(i = 0; i < 2; i++) {
181 LocDim[i] = 0.0;
182 LocIx[i] = 0.0;
183 LocIy[i] = 0.0;
184 LocIz[i] = 0.0;
185 LocIxx[i] = 0.0;
186 LocIyy[i] = 0.0;
187 LocIzz[i] = 0.0;
188 LocIxy[i] = 0.0;
189 LocIyz[i] = 0.0;
190 LocIxz[i] = 0.0;
191 }
192
193 NbUGaussP[0] = S.SIntOrder(EpsDim);
194 NbUGaussP[1] = RealToInt(Ceiling(ERROR_ALGEBR_RATIO*NbUGaussP[0]));
195 math::GaussPoints(NbUGaussP[0],UGaussP0); math::GaussWeights(NbUGaussP[0],UGaussW0);
196 math::GaussPoints(NbUGaussP[1],UGaussP1); math::GaussWeights(NbUGaussP[1],UGaussW1);
197
198 NbUSubs = S.SUIntSubs();
199 TColStd_Array1OfReal UKnots(1,NbUSubs+1);
200 S.UKnots(UKnots);
201
202 while (isNaturalRestriction || D.More()) {
203 if(isNaturalRestriction){
204 NbLGaussP[0] = Min(2*NbUGaussP[0],math::GaussPointsMax());
205 }else{
206 S.Load(D.Value()); ++iD;
207 NbLGaussP[0] = S.LIntOrder(EpsDim);
208 }
209 NbLGaussP[1] = RealToInt(Ceiling(ERROR_ALGEBR_RATIO*NbLGaussP[0]));
210 math::GaussPoints(NbLGaussP[0],LGaussP0); math::GaussWeights(NbLGaussP[0],LGaussW0);
211 math::GaussPoints(NbLGaussP[1],LGaussP1); math::GaussWeights(NbLGaussP[1],LGaussW1);
212
213 NbLSubs = isNaturalRestriction? S.SVIntSubs(): S.LIntSubs();
214 TColStd_Array1OfReal LKnots(1,NbLSubs+1);
215 if(isNaturalRestriction){
216 S.VKnots(LKnots);
217 l1 = BV1; l2 = BV2;
218 }else{
219 S.LKnots(LKnots);
220 l1 = S.FirstParameter(); l2 = S.LastParameter();
221 }
222 ErrorL = 0.0;
223 kLEnd = 1; JL = 0;
224 //OCC503(apo): if(Abs(l2-l1) < EPS_PARAM) continue;
225 if(Abs(l2-l1) > EPS_PARAM) {
226 iLSubEnd = LFillIntervalBounds(l1, l2, LKnots, NumSubs);
227 LMaxSubs = MaxSubs(iLSubEnd);
228 //-- exception avoiding
229 if(LMaxSubs > SM) LMaxSubs = SM;
230 DimL.Init(0.0,1,LMaxSubs); ErrL.Init(0.0,1,LMaxSubs); ErrUL.Init(0.0,1,LMaxSubs);
231 do{// while: L
232 if(++JL > iLSubEnd){
233 LRange[0] = IL = ErrL->Max(); LRange[1] = JL;
234 L1(JL) = (L1(IL) + L2(IL))/2.0; L2(JL) = L2(IL); L2(IL) = L1(JL);
235 }else LRange[0] = IL = JL;
236 if(JL == LMaxSubs || Abs(L2(JL) - L1(JL)) < EPS_PARAM)
237 if(kLEnd == 1){
238 DimL(JL) = ErrL(JL) = IxL(JL) = IyL(JL) = IzL(JL) =
239 IxxL(JL) = IyyL(JL) = IzzL(JL) = IxyL(JL) = IxzL(JL) = IyzL(JL) = 0.0;
240 }else{
241 JL--;
242 EpsL = ErrorL; Eps = EpsL/0.9;
243 break;
244 }
245 else
246 for(kL=0; kL < kLEnd; kL++){
247 iLS = LRange[kL];
248 lm = 0.5*(L2(iLS) + L1(iLS));
249 lr = 0.5*(L2(iLS) - L1(iLS));
250 CIx = CIy = CIz = CIxy = CIxz = CIyz = 0.0;
251 for(iGL=0; iGL < iGLEnd; iGL++){//
252 CDim[iGL] = CIxx[iGL] = CIyy[iGL] = CIzz[iGL] = 0.0;
253 for(iL=1; iL<=NbLGaussP[iGL]; iL++){
254 l = lm + lr*(*LGaussP[iGL])(iL);
255 if(isNaturalRestriction){
256 v = l; u2 = BU2; Dul = (*LGaussW[iGL])(iL);
257 }else{
258 S.D12d (l, Puv, Vuv);
259 Dul = Vuv.Y()*(*LGaussW[iGL])(iL); // Dul = Du / Dl
260 if(Abs(Dul) < EPS_PARAM) continue;
261 v = Puv.Y(); u2 = Puv.X();
262 //Check on cause out off bounds of value current parameter
263 if(v < BV1) v = BV1; else if(v > BV2) v = BV2;
264 if(u2 < BU1) u2 = BU1; else if(u2 > BU2) u2 = BU2;
265 }
266 ErrUL(iLS) = 0.0;
267 kUEnd = 1; JU = 0;
268 if(Abs(u2-u1) < EPS_PARAM) continue;
269 iUSubEnd = UFillIntervalBounds(u1, u2, UKnots, NumSubs);
270 UMaxSubs = MaxSubs(iUSubEnd);
271 //-- exception avoiding
272 if(UMaxSubs > SM) UMaxSubs = SM;
273 DimU.Init(0.0,1,UMaxSubs); ErrU.Init(0.0,1,UMaxSubs); ErrorU = 0.0;
274 do{//while: U
275 if(++JU > iUSubEnd){
276 URange[0] = IU = ErrU->Max(); URange[1] = JU;
277 U1(JU) = (U1(IU)+U2(IU))/2.0; U2(JU) = U2(IU); U2(IU) = U1(JU);
278 }else URange[0] = IU = JU;
279 if(JU == UMaxSubs || Abs(U2(JU) - U1(JU)) < EPS_PARAM)
280 if(kUEnd == 1){
281 DimU(JU) = ErrU(JU) = IxU(JU) = IyU(JU) = IzU(JU) =
282 IxxU(JU) = IyyU(JU) = IzzU(JU) = IxyU(JU) = IxzU(JU) = IyzU(JU) = 0.0;
283 }else{
284 JU--;
285 EpsU = ErrorU; Eps = EpsU*Abs((u2-u1)*Dul)/0.1; EpsL = 0.9*Eps;
286 break;
287 }
288 else
289 for(kU=0; kU < kUEnd; kU++){
290 iUS = URange[kU];
291 um = 0.5*(U2(iUS) + U1(iUS));
292 ur = 0.5*(U2(iUS) - U1(iUS));
293 iGUEnd = iGLEnd - iGL;
294 for(iGU=0; iGU < iGUEnd; iGU++){//
295 LocDim[iGU] =
296 LocIxx[iGU] = LocIyy[iGU] = LocIzz[iGU] =
297 LocIx[iGU] = LocIy[iGU] = LocIz[iGU] =
298 LocIxy[iGU] = LocIxz[iGU] = LocIyz[iGU] = 0.0;
299 for(iU=1; iU<=NbUGaussP[iGU]; iU++){
300 u = um + ur*(*UGaussP[iGU])(iU);
301 S.Normal(u, v, Ps, VNor);
302 VNor.Coord(xn, yn, zn);
303 Ps.Coord(x, y, z);
304 x -= xloc; y -= yloc; z -= zloc;
305 xn *= (*UGaussW[iGU])(iU);
306 yn *= (*UGaussW[iGU])(iU);
307 zn *= (*UGaussW[iGU])(iU);
308 if(ByPoint){
309 //volume of elementary cone
310 dDim = (x*xn+y*yn+z*zn)/3.0;
311 //coordinates of cone's center mass
312 px = 0.75*x; py = 0.75*y; pz = 0.75*z;
313 LocDim[iGU] += dDim;
314 //if(iGU > 0) continue;
315 LocIx[iGU] += px*dDim;
316 LocIy[iGU] += py*dDim;
317 LocIz[iGU] += pz*dDim;
318 x -= Coeff[0]; y -= Coeff[1]; z -= Coeff[2];
319 dDim *= 3.0/5.0;
320 LocIxy[iGU] -= x*y*dDim;
321 LocIyz[iGU] -= y*z*dDim;
322 LocIxz[iGU] -= x*z*dDim;
323 xi = x*x; yi = y*y; zi = z*z;
324 LocIxx[iGU] += (yi+zi)*dDim;
325 LocIyy[iGU] += (xi+zi)*dDim;
326 LocIzz[iGU] += (xi+yi)*dDim;
327 }else{ // by plane
328 s = xn*Coeff[0] + yn*Coeff[1] + zn*Coeff[2];
329 d1 = Coeff[0]*x + Coeff[1]*y + Coeff[2]*z - Coeff[3];
330 d2 = d1*d1;
331 d3 = d1*d2/3.0;
332 ds = s*d1;
333 LocDim[iGU] += ds;
334 //if(iGU > 0) continue;
335 LocIx[iGU] += (x - Coeff[0]*d1/2.0) * ds;
336 LocIy[iGU] += (y - Coeff[1]*d1/2.0) * ds;
337 LocIz[iGU] += (z - Coeff[2]*d1/2.0) * ds;
338 px = x-Coeff[0]*d1; py = y-Coeff[1]*d1; pz = z-Coeff[2]*d1;
339 xi = px*px*d1 + px*Coeff[0]*d2 + Coeff[0]*Coeff[0]*d3;
340 yi = py*py*d1 + py*Coeff[1]*d2 + Coeff[1]*Coeff[1]*d3;
341 zi = pz*pz*d1 + pz*Coeff[2]*d2 + Coeff[2]*Coeff[2]*d3;
342 LocIxx[iGU] += (yi+zi)*s;
343 LocIyy[iGU] += (xi+zi)*s;
344 LocIzz[iGU] += (xi+yi)*s;
345 d2 /= 2.0;
346 xi = py*pz*d1 + py*Coeff[2]*d2 + pz*Coeff[1]*d2 + Coeff[1]*Coeff[2]*d3;
347 yi = px*pz*d1 + pz*Coeff[0]*d2 + px*Coeff[2]*d2 + Coeff[0]*Coeff[2]*d3;
348 zi = px*py*d1 + px*Coeff[1]*d2 + py*Coeff[0]*d2 + Coeff[0]*Coeff[1]*d3;
349 LocIxy[iGU] -= zi*s; LocIyz[iGU] -= xi*s; LocIxz[iGU] -= yi*s;
350 }
351 }//for: iU
352 }//for: iGU
353 DimU(iUS) = LocDim[0]*ur;
354 IxxU(iUS) = LocIxx[0]*ur; IyyU(iUS) = LocIyy[0]*ur; IzzU(iUS) = LocIzz[0]*ur;
355 if(iGL > 0) continue;
356 LocDim[1] = Abs(LocDim[1]-LocDim[0]);
357 LocIxx[1] = Abs(LocIxx[1]-LocIxx[0]);
358 LocIyy[1] = Abs(LocIyy[1]-LocIyy[0]);
359 LocIzz[1] = Abs(LocIzz[1]-LocIzz[0]);
360 ErrU(iUS) = isMinDim? LocDim[1]*ur: (LocIxx[1] + LocIyy[1] + LocIzz[1])*ur;
361 IxU(iUS) = LocIx[0]*ur; IyU(iUS) = LocIy[0]*ur; IzU(iUS) = LocIz[0]*ur;
362 IxyU(iUS) = LocIxy[0]*ur; IxzU(iUS) = LocIxz[0]*ur; IyzU(iUS) = LocIyz[0]*ur;
363 }//for: kU (iUS)
364 if(JU == iUSubEnd) kUEnd = 2;
365 if(kUEnd == 2) {
366 Standard_Integer imax = ErrU->Max();
367 if(imax > 0) ErrorU = ErrU(imax);
368 else ErrorU = 0.0;
369 }
370 }while((ErrorU - EpsU > 0.0 && EpsU != 0.0) || kUEnd == 1);
371 for(i=1; i<=JU; i++) {
372 CDim[iGL] += DimU(i)*Dul;
373 CIxx[iGL] += IxxU(i)*Dul; CIyy[iGL] += IyyU(i)*Dul; CIzz[iGL] += IzzU(i)*Dul;
374 }
375 if(iGL > 0) continue;
376 ErrUL(iLS) = ErrorU*Abs((u2-u1)*Dul);
377 for(i=1; i<=JU; i++){
378 CIx += IxU(i)*Dul; CIy += IyU(i)*Dul; CIz += IzU(i)*Dul;
379 //CIxx += IxxU(i)*Dul; CIyy += IyyU(i)*Dul; CIzz += IzzU(i)*Dul;
380 CIxy += IxyU(i)*Dul; CIxz += IxzU(i)*Dul; CIyz += IyzU(i)*Dul;
381 }
382 }//for: iL
383 }//for: iGL
384 DimL(iLS) = CDim[0]*lr;
385 IxxL(iLS) = CIxx[0]*lr; IyyL(iLS) = CIyy[0]*lr; IzzL(iLS) = CIzz[0]*lr;
386 if(iGLEnd == 2) {
387 //ErrL(iLS) = Abs(CDim[1]-CDim[0])*lr + ErrUL(iLS);
388 CDim[1] = Abs(CDim[1]-CDim[0]);
389 CIxx[1] = Abs(CIxx[1]-CIxx[0]); CIyy[1] = Abs(CIyy[1]-CIyy[0]); CIzz[1] = Abs(CIzz[1]-CIzz[0]);
390 ErrorU = ErrUL(iLS);
391 ErrL(iLS) = (isMinDim? CDim[1]: (CIxx[1] + CIyy[1] + CIzz[1]))*lr + ErrorU;
392 }
393 IxL(iLS) = CIx*lr; IyL(iLS) = CIy*lr; IzL(iLS) = CIz*lr;
394 //IxxL(iLS) = CIxx*lr; IyyL(iLS) = CIyy*lr; IzzL(iLS) = CIzz*lr;
395 IxyL(iLS) = CIxy*lr; IxzL(iLS) = CIxz*lr; IyzL(iLS) = CIyz*lr;
396 }//for: (kL)iLS
397 // Calculate/correct epsilon of computation by current value of Dim
398 //That is need for not spend time for
399 if(JL == iLSubEnd){
400 kLEnd = 2;
401 Standard_Real DDim = 0.0, DIxx = 0.0, DIyy = 0.0, DIzz = 0.0;
402 for(i=1; i<=JL; i++) {
403 DDim += DimL(i);
404 DIxx += IxxL(i); DIyy += IyyL(i); DIzz += IzzL(i);
405 }
406 DDim = isMinDim? Abs(DDim): Abs(DIxx) + Abs(DIyy) + Abs(DIzz);
407 DDim = Abs(DDim*EpsDim);
408 if(DDim > Eps) {
409 Eps = DDim; EpsL = 0.9*Eps;
410 }
411 }
412 if(kLEnd == 2) {
413 Standard_Integer imax = ErrL->Max();
414 if(imax > 0) ErrorL = ErrL(imax);
415 else ErrorL = 0.0;
416 }
417 }while((ErrorL - EpsL > 0.0 && isVerifyComputation) || kLEnd == 1);
418 for(i=1; i<=JL; i++){
419 Dim += DimL(i);
420 Ix += IxL(i); Iy += IyL(i); Iz += IzL(i);
421 Ixx += IxxL(i); Iyy += IyyL(i); Izz += IzzL(i);
422 Ixy += IxyL(i); Ixz += IxzL(i); Iyz += IyzL(i);
423 }
424 ErrorLMax = Max(ErrorLMax, ErrorL);
425 }
426 if(isNaturalRestriction) break;
427 D.Next();
428 }
429 if(Abs(Dim) >= EPS_DIM){
430 if(ByPoint){
431 Ix = Coeff[0] + Ix/Dim;
432 Iy = Coeff[1] + Iy/Dim;
433 Iz = Coeff[2] + Iz/Dim;
434 }else{
435 Ix /= Dim;
436 Iy /= Dim;
437 Iz /= Dim;
438 }
439 g.SetCoord (Ix, Iy, Iz);
440 }else{
441 Dim =0.;
442 g.SetCoord(0.,0.,0.);
443 }
444 inertia.SetCols (gp_XYZ (Ixx, Ixy, Ixz),
445 gp_XYZ (Ixy, Iyy, Iyz),
446 gp_XYZ (Ixz, Iyz, Izz));
447 if(iGLEnd == 2)
448 Eps = Dim != 0.0? ErrorLMax/(isMinDim? Abs(Dim): (Abs(Ixx) + Abs(Iyy) + Abs(Izz))): 0.0;
449 else Eps = EpsDim;
450 return Eps;
451}
452
453static Standard_Real Compute(Face& S, const Standard_Boolean ByPoint, const Standard_Real Coeff[],
454 const gp_Pnt& loc, Standard_Real& Dim, gp_Pnt& g, gp_Mat& inertia, Standard_Real EpsDim)
455{
456 Standard_Boolean isErrorCalculation = 0.0 > EpsDim || EpsDim < 0.001? 1: 0;
457 Standard_Boolean isVerifyComputation = 0.0 < EpsDim && EpsDim < 0.001? 1: 0;
458 EpsDim = Abs(EpsDim);
459 Domain D;
460 return CCompute(S,D,ByPoint,Coeff,loc,Dim,g,inertia,EpsDim,isErrorCalculation,isVerifyComputation);
461}
462
463static Standard_Real Compute(Face& S, Domain& D, const Standard_Boolean ByPoint, const Standard_Real Coeff[],
464 const gp_Pnt& loc, Standard_Real& Dim, gp_Pnt& g, gp_Mat& inertia, Standard_Real EpsDim)
465{
466 Standard_Boolean isErrorCalculation = 0.0 > EpsDim || EpsDim < 0.001? 1: 0;
467 Standard_Boolean isVerifyComputation = 0.0 < EpsDim && EpsDim < 0.001? 1: 0;
468 EpsDim = Abs(EpsDim);
469 return CCompute(S,D,ByPoint,Coeff,loc,Dim,g,inertia,EpsDim,isErrorCalculation,isVerifyComputation);
470}
471
472static void Compute(const Face& S,
473 const Standard_Boolean ByPoint,
474 const Standard_Real Coeff[],
475 const gp_Pnt& Loc,
476 Standard_Real& Volu,
477 gp_Pnt& G,
478 gp_Mat& Inertia)
479{
480
481 gp_Pnt P;
482 gp_Vec VNor;
483 Standard_Real dvi, dv;
484 Standard_Real ur, um, u, vr, vm, v;
485 Standard_Real x, y, z, xn, yn, zn, xi, yi, zi;
486// Standard_Real x, y, z, xn, yn, zn, xi, yi, zi, xyz;
487 Standard_Real px,py,pz,s,d1,d2,d3;
488 Standard_Real Ixi, Iyi, Izi, Ixxi, Iyyi, Izzi, Ixyi, Ixzi, Iyzi;
489 Standard_Real xloc, yloc, zloc;
490 Standard_Real Ix, Iy, Iz, Ixx, Iyy, Izz, Ixy, Ixz, Iyz;
491
492 Volu = Ix = Iy = Iz = Ixx = Iyy = Izz = Ixy = Ixz = Iyz = 0.0;
493 Loc.Coord (xloc, yloc, zloc);
494
495 Standard_Real LowerU, UpperU, LowerV, UpperV;
496 S.Bounds ( LowerU, UpperU, LowerV, UpperV);
497 Standard_Integer UOrder = Min(S.UIntegrationOrder (),
498 math::GaussPointsMax());
499 Standard_Integer VOrder = Min(S.VIntegrationOrder (),
500 math::GaussPointsMax());
501
502 Standard_Integer i, j;
503 math_Vector GaussPU (1, UOrder); //gauss points and weights
504 math_Vector GaussWU (1, UOrder);
505 math_Vector GaussPV (1, VOrder);
506 math_Vector GaussWV (1, VOrder);
507
508 math::GaussPoints (UOrder,GaussPU);
509 math::GaussWeights (UOrder,GaussWU);
510 math::GaussPoints (VOrder,GaussPV);
511 math::GaussWeights (VOrder,GaussWV);
512
513 um = 0.5 * (UpperU + LowerU);
514 vm = 0.5 * (UpperV + LowerV);
515 ur = 0.5 * (UpperU - LowerU);
516 vr = 0.5 * (UpperV - LowerV);
517
518 for (j = 1; j <= VOrder; j++) {
519 v = vm + vr * GaussPV (j);
520 dvi = Ixi = Iyi = Izi = Ixxi = Iyyi = Izzi = Ixyi = Ixzi = Iyzi = 0.0;
521
522 for (i = 1; i <= UOrder; i++) {
523 u = um + ur * GaussPU (i);
524 S.Normal (u, v, P, VNor);
525 VNor.Coord (xn, yn, zn);
526 P.Coord (x, y, z);
527 x -= xloc; y -= yloc; z -= zloc;
528 xn *= GaussWU (i); yn *= GaussWU (i); zn *= GaussWU (i);
529 if (ByPoint) {
530 ///////////////////// ///////////////////////
531 // OFV code // // Initial code //
532 ///////////////////// ///////////////////////
533 // modified by APO
534 dv = (x*xn+y*yn+z*zn)/3.0; //xyz = x * y * z;
535 dvi += dv; //Ixyi += zn * xyz;
536 Ixi += 0.75*x*dv; //Iyzi += xn * xyz;
537 Iyi += 0.75*y*dv; //Ixzi += yn * xyz;
538 Izi += 0.75*z*dv; //xi = x * x * x * xn / 3.0;
539 x -= Coeff[0]; //yi = y * y * y * yn / 3.0;
540 y -= Coeff[1]; //zi = z * z * z * zn / 3.0;
541 z -= Coeff[2]; //Ixxi += (yi + zi);
542 dv *= 3.0/5.0; //Iyyi += (xi + zi);
543 Ixyi -= x*y*dv; //Izzi += (xi + yi);
544 Iyzi -= y*z*dv; //x -= Coeff[0];
545 Ixzi -= x*z*dv; //y -= Coeff[1];
546 xi = x*x; //z -= Coeff[2];
547 yi = y*y; //dv = x * xn + y * yn + z * zn;
548 zi = z*z; //dvi += dv;
549 Ixxi += (yi + zi)*dv; //Ixi += x * dv;
550 Iyyi += (xi + zi)*dv; //Iyi += y * dv;
551 Izzi += (xi + yi)*dv; //Izi += z * dv;
552 }
553 else { // by plane
554 s = xn * Coeff[0] + yn * Coeff[1] + zn * Coeff[2];
555 d1 = Coeff[0] * x + Coeff[1] * y + Coeff[2] * z - Coeff[3];
556 d2 = d1 * d1;
557 d3 = d1 * d2 / 3.0;
558 dv = s * d1;
559 dvi += dv;
560 Ixi += (x - (Coeff[0] * d1 / 2.0)) * dv;
561 Iyi += (y - (Coeff[1] * d1 / 2.0)) * dv;
562 Izi += (z - (Coeff[2] * d1 / 2.0)) * dv;
563 px = x - Coeff[0] * d1;
564 py = y - Coeff[1] * d1;
565 pz = z - Coeff[2] * d1;
566 xi = px * px * d1 + px * Coeff[0]* d2 + Coeff[0] * Coeff[0] * d3;
567 yi = py * py * d1 + py * Coeff[1] * d2 + Coeff[1] * Coeff[1] * d3;
568 zi = pz * pz * d1 + pz * Coeff[2] * d2 + Coeff[2] * Coeff[2] * d3;
569 Ixxi += (yi + zi) * s;
570 Iyyi += (xi + zi) * s;
571 Izzi += (xi + yi) * s;
572 d2 /= 2.0;
573 xi = (py * pz * d1) + (py * Coeff[2] * d2) + (pz * Coeff[1] * d2) + (Coeff[1] * Coeff[2] * d3);
574 yi = (px * pz * d1) + (pz * Coeff[0] * d2) + (px * Coeff[2] * d2) + (Coeff[0] * Coeff[2] * d3);
575 zi = (px * py * d1) + (px * Coeff[1] * d2) + (py * Coeff[0] * d2) + (Coeff[0] * Coeff[1] * d3);
576 Ixyi -= zi * s;
577 Iyzi -= xi * s;
578 Ixzi -= yi * s;
579 }
580 }
581 Volu += dvi * GaussWV (j);
582 Ix += Ixi * GaussWV (j);
583 Iy += Iyi * GaussWV (j);
584 Iz += Izi * GaussWV (j);
585 Ixx += Ixxi * GaussWV (j);
586 Iyy += Iyyi * GaussWV (j);
587 Izz += Izzi * GaussWV (j);
588 Ixy += Ixyi * GaussWV (j);
589 Ixz += Ixzi * GaussWV (j);
590 Iyz += Iyzi * GaussWV (j);
591 }
592 vr *= ur;
593 Ixx *= vr;
594 Iyy *= vr;
595 Izz *= vr;
596 Ixy *= vr;
597 Ixz *= vr;
598 Iyz *= vr;
599 if (Abs(Volu) >= EPS_DIM ) {
600 if (ByPoint) {
601 Ix = Coeff[0] + Ix/Volu;
602 Iy = Coeff[1] + Iy/Volu;
603 Iz = Coeff[2] + Iz/Volu;
604 Volu *= vr;
605 }
606 else { //by plane
607 Ix /= Volu;
608 Iy /= Volu;
609 Iz /= Volu;
610 Volu *= vr;
611 }
612 G.SetCoord (Ix, Iy, Iz);
613 }
614 else {
615 G.SetCoord(0.,0.,0.);
616 Volu =0.;
617 }
618 Inertia.SetCols (gp_XYZ (Ixx, Ixy, Ixz),
619 gp_XYZ (Ixy, Iyy, Iyz),
620 gp_XYZ (Ixz, Iyz, Izz));
621
622}
623
624// Last modified by OFV 5.2001:
625// 1). surface and edge integration order is equal now
626// 2). "by point" works now rathre correctly (it looks so...)
627static void Compute(Face& S, Domain& D, const Standard_Boolean ByPoint, const Standard_Real Coeff[],
628 const gp_Pnt& Loc, Standard_Real& Volu, gp_Pnt& G, gp_Mat& Inertia)
629
630{
631 Standard_Real x, y, z, xi, yi, zi, l1, l2, lm, lr, l, v1, v2, v, u1, u2, um, ur, u, ds, Dul, xloc, yloc, zloc;
632 Standard_Real LocVolu, LocIx, LocIy, LocIz, LocIxx, LocIyy, LocIzz, LocIxy, LocIxz, LocIyz;
633 Standard_Real CVolu, CIx, CIy, CIz, CIxx, CIyy, CIzz, CIxy, CIxz, CIyz, Ix, Iy, Iz, Ixx, Iyy, Izz, Ixy, Ixz, Iyz;
634 Standard_Real xn, yn, zn, px, py, pz, s, d1, d2, d3, dSigma;
635 Standard_Integer i, j, vio, sio, max, NbGaussgp_Pnts;
636
637 gp_Pnt Ps;
638 gp_Vec VNor;
639 gp_Pnt2d Puv;
640 gp_Vec2d Vuv;
641
642 Loc.Coord (xloc, yloc, zloc);
643 Volu = Ix = Iy = Iz = Ixx = Iyy = Izz = Ixy = Ixz = Iyz = 0.0;
644 S.Bounds (u1, u2, v1, v2);
645 Standard_Real _u2 = u2; //OCC104
646 vio = S.VIntegrationOrder ();
647
648 while (D.More())
649 {
650 S.Load(D.Value());
651 sio = S.IntegrationOrder ();
652 max = Max(vio,sio);
653 NbGaussgp_Pnts = Min(max,math::GaussPointsMax());
654
655 math_Vector GaussP (1, NbGaussgp_Pnts);
656 math_Vector GaussW (1, NbGaussgp_Pnts);
657 math::GaussPoints (NbGaussgp_Pnts,GaussP);
658 math::GaussWeights (NbGaussgp_Pnts,GaussW);
659
660 CVolu = CIx = CIy = CIz = CIxx = CIyy = CIzz = CIxy = CIxz = CIyz = 0.0;
661 l1 = S.FirstParameter();
662 l2 = S.LastParameter();
663 lm = 0.5 * (l2 + l1);
664 lr = 0.5 * (l2 - l1);
665
666 for (i=1; i<=NbGaussgp_Pnts; i++)
667 {
668 l = lm + lr * GaussP(i);
669 S.D12d (l, Puv, Vuv);
670 v = Puv.Y();
671 u2 = Puv.X();
672
673 //OCC104
674 v = v < v1? v1: v;
675 v = v > v2? v2: v;
676 u2 = u2 < u1? u1: u2;
677 u2 = u2 > _u2? _u2: u2;
678
679 Dul = Vuv.Y() * GaussW(i);
680 um = 0.5 * (u2 + u1);
681 ur = 0.5 * (u2 - u1);
682 LocVolu = LocIx = LocIy = LocIz = LocIxx = LocIyy = LocIzz = LocIxy = LocIxz = LocIyz = 0.0;
683
684 for (j=1; j<=NbGaussgp_Pnts; j++)
685 {
686 u = um + ur * GaussP(j);
687 S.Normal (u, v, Ps, VNor);
688 VNor.Coord (xn, yn, zn);
689 Ps.Coord (x, y, z);
690 x -= xloc;
691 y -= yloc;
692 z -= zloc;
693 xn = xn * Dul * GaussW(j);
694 yn = yn * Dul * GaussW(j);
695 zn = zn * Dul * GaussW(j);
696 if(ByPoint)
697 {
698 dSigma = (x*xn+y*yn+z*zn)/3.0;
699 LocVolu += dSigma;
700 LocIx += 0.75*x*dSigma;
701 LocIy += 0.75*y*dSigma;
702 LocIz += 0.75*z*dSigma;
703 x -= Coeff[0];
704 y -= Coeff[1];
705 z -= Coeff[2];
706 dSigma *= 3.0/5.0;
707 LocIxy -= x*y*dSigma;
708 LocIyz -= y*z*dSigma;
709 LocIxz -= x*z*dSigma;
710 xi = x*x;
711 yi = y*y;
712 zi = z*z;
713 LocIxx += (yi + zi)*dSigma;
714 LocIyy += (xi + zi)*dSigma;
715 LocIzz += (xi + yi)*dSigma;
716 }
717 else
718 {
719 s = xn * Coeff[0] + yn * Coeff[1] + zn * Coeff[2];
720 d1 = Coeff[0] * x + Coeff[1] * y + Coeff[2] * z;
721 d2 = d1 * d1;
722 d3 = d1 * d2 / 3.0;
723 ds = s * d1;
724 LocVolu += ds;
725 LocIx += (x - Coeff[0] * d1 / 2.0) * ds;
726 LocIy += (y - Coeff[1] * d1 / 2.0) * ds;
727 LocIz += (z - Coeff[2] * d1 / 2.0) * ds;
728 px = x - Coeff[0] * d1;
729 py = y - Coeff[1] * d1;
730 pz = z - Coeff[2] * d1;
731 xi = (px * px * d1) + (px * Coeff[0]* d2) + (Coeff[0] * Coeff[0] * d3);
732 yi = (py * py * d1) + (py * Coeff[1] * d2) + (Coeff[1] * Coeff[1] * d3);
733 zi = pz * pz * d1 + pz * Coeff[2] * d2 + (Coeff[2] * Coeff[2] * d3);
734 LocIxx += (yi + zi) * s;
735 LocIyy += (xi + zi) * s;
736 LocIzz += (xi + yi) * s;
737 d2 /= 2.0;
738 xi = (py * pz * d1) + (py * Coeff[2] * d2) + (pz * Coeff[1] * d2) + (Coeff[1] * Coeff[2] * d3);
739 yi = (px * pz * d1) + (pz * Coeff[0] * d2) + (px * Coeff[2] * d2) + (Coeff[0] * Coeff[2] * d3);
740 zi = (px * py * d1) + (px * Coeff[1] * d2) + (py * Coeff[0] * d2) + (Coeff[0] * Coeff[1] * d3);
741 LocIxy -= zi * s;
742 LocIyz -= xi * s;
743 LocIxz -= yi * s;
744 }
745 }
746 CVolu += LocVolu * ur;
747 CIx += LocIx * ur;
748 CIy += LocIy * ur;
749 CIz += LocIz * ur;
750 CIxx += LocIxx * ur;
751 CIyy += LocIyy * ur;
752 CIzz += LocIzz * ur;
753 CIxy += LocIxy * ur;
754 CIxz += LocIxz * ur;
755 CIyz += LocIyz * ur;
756 }
757 Volu += CVolu * lr;
758 Ix += CIx * lr;
759 Iy += CIy * lr;
760 Iz += CIz * lr;
761 Ixx += CIxx * lr;
762 Iyy += CIyy * lr;
763 Izz += CIzz * lr;
764 Ixy += CIxy * lr;
765 Ixz += CIxz * lr;
766 Iyz += CIyz * lr;
767 D.Next();
768 }
769
770 if(Abs(Volu) >= EPS_DIM)
771 {
772 if(ByPoint)
773 {
774 Ix = Coeff[0] + Ix/Volu;
775 Iy = Coeff[1] + Iy/Volu;
776 Iz = Coeff[2] + Iz/Volu;
777 }
778 else
779 {
780 Ix /= Volu;
781 Iy /= Volu;
782 Iz /= Volu;
783 }
784 G.SetCoord (Ix, Iy, Iz);
785 }
786 else
787 {
788 Volu =0.;
789 G.SetCoord(0.,0.,0.);
790 }
791
792 Inertia.SetCols (gp_XYZ (Ixx, Ixy, Ixz),
793 gp_XYZ (Ixy, Iyy, Iyz),
794 gp_XYZ (Ixz, Iyz, Izz));
795
796}
797
798GProp_VGProps::GProp_VGProps(){}
799
800GProp_VGProps::GProp_VGProps(Face& S, const gp_Pnt& VLocation, const Standard_Real Eps){
801 SetLocation(VLocation);
802 Perform(S,Eps);
803}
804
805GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pnt& VLocation, const Standard_Real Eps){
806 SetLocation(VLocation);
807 Perform(S,D,Eps);
808}
809
810GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pnt& VLocation){
811 SetLocation(VLocation);
812 Perform(S,D);
813}
814
815GProp_VGProps::GProp_VGProps(const Face& S, const gp_Pnt& VLocation){
816 SetLocation(VLocation);
817 Perform(S);
818}
819
820GProp_VGProps::GProp_VGProps(Face& S, const gp_Pnt& O, const gp_Pnt& VLocation, const Standard_Real Eps){
821 SetLocation(VLocation);
822 Perform(S,O,Eps);
823}
824
825GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pnt& O, const gp_Pnt& VLocation, const Standard_Real Eps){
826 SetLocation(VLocation);
827 Perform(S,D,O,Eps);
828}
829
830GProp_VGProps::GProp_VGProps(const Face& S, const gp_Pnt& O, const gp_Pnt& VLocation){
831 SetLocation(VLocation);
832 Perform(S,O);
833}
834
835GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pnt& O, const gp_Pnt& VLocation){
836 SetLocation(VLocation);
837 Perform(S,D,O);
838}
839
840GProp_VGProps::GProp_VGProps(Face& S, const gp_Pln& Pl, const gp_Pnt& VLocation, const Standard_Real Eps){
841 SetLocation(VLocation);
842 Perform(S,Pl,Eps);
843}
844
845GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pln& Pl, const gp_Pnt& VLocation, const Standard_Real Eps){
846 SetLocation(VLocation);
847 Perform(S,D,Pl,Eps);
848}
849
850GProp_VGProps::GProp_VGProps(const Face& S, const gp_Pln& Pl, const gp_Pnt& VLocation){
851 SetLocation(VLocation);
852 Perform(S,Pl);
853}
854
855GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pln& Pl, const gp_Pnt& VLocation){
856 SetLocation(VLocation);
857 Perform(S,D,Pl);
858}
859
860void GProp_VGProps::SetLocation(const gp_Pnt& VLocation){
861 loc = VLocation;
862}
863
864Standard_Real GProp_VGProps::Perform(Face& S, const Standard_Real Eps){
865 Standard_Real Coeff[] = {0., 0., 0.};
866 return myEpsilon = Compute(S,Standard_True,Coeff,loc,dim,g,inertia,Eps);
867}
868
869Standard_Real GProp_VGProps::Perform(Face& S, Domain& D, const Standard_Real Eps){
870 Standard_Real Coeff[] = {0., 0., 0.};
871 return myEpsilon = Compute(S,D,Standard_True,Coeff,loc,dim,g,inertia,Eps);
872}
873
874void GProp_VGProps::Perform(const Face& S){
875 Standard_Real Coeff[] = {0., 0., 0.};
876 Compute(S,Standard_True,Coeff,loc,dim,g,inertia);
877 myEpsilon = 1.0;
878 return;
879}
880
881void GProp_VGProps::Perform(Face& S, Domain& D){
882 Standard_Real Coeff[] = {0., 0., 0.};
883 Compute(S,D,Standard_True,Coeff,loc,dim,g,inertia);
884 myEpsilon = 1.0;
885 return;
886}
887
888Standard_Real GProp_VGProps::Perform(Face& S, const gp_Pnt& O, const Standard_Real Eps){
889 Standard_Real xloc, yloc, zloc;
890 loc.Coord(xloc, yloc, zloc);
891 Standard_Real Coeff[3];
892 O.Coord (Coeff[0], Coeff[1], Coeff[2]);
893 Coeff[0] -= xloc; Coeff[1] -= yloc; Coeff[2] -= zloc;
894 return myEpsilon = Compute(S,Standard_True,Coeff,loc,dim,g,inertia,Eps);
895}
896
897Standard_Real GProp_VGProps::Perform(Face& S, Domain& D, const gp_Pnt& O, const Standard_Real Eps){
898 Standard_Real xloc, yloc, zloc;
899 loc.Coord(xloc, yloc, zloc);
900 Standard_Real Coeff[3];
901 O.Coord (Coeff[0], Coeff[1], Coeff[2]);
902 Coeff[0] -= xloc; Coeff[1] -= yloc; Coeff[2] -= zloc;
903 return myEpsilon = Compute(S,D,Standard_True,Coeff,loc,dim,g,inertia,Eps);
904}
905
906void GProp_VGProps::Perform(const Face& S, const gp_Pnt& O){
907 Standard_Real xloc, yloc, zloc;
908 loc.Coord(xloc, yloc, zloc);
909 Standard_Real Coeff[3];
910 O.Coord (Coeff[0], Coeff[1], Coeff[2]);
911 Coeff[0] -= xloc; Coeff[1] -= yloc; Coeff[2] -= zloc;
912 Compute(S,Standard_True,Coeff,loc,dim,g,inertia);
913 myEpsilon = 1.0;
914 return;
915}
916
917void GProp_VGProps::Perform(Face& S, Domain& D, const gp_Pnt& O){
918 Standard_Real xloc, yloc, zloc;
919 loc.Coord(xloc, yloc, zloc);
920 Standard_Real Coeff[3];
921 O.Coord (Coeff[0], Coeff[1], Coeff[2]);
922 Coeff[0] -= xloc; Coeff[1] -= yloc; Coeff[2] -= zloc;
923 Compute(S,D,Standard_True,Coeff,loc,dim,g,inertia);
924 myEpsilon = 1.0;
925 return;
926}
927
928Standard_Real GProp_VGProps::Perform(Face& S, const gp_Pln& Pl, const Standard_Real Eps){
929 Standard_Real xloc, yloc, zloc;
930 loc.Coord (xloc, yloc, zloc);
931 Standard_Real Coeff[4];
932 Pl.Coefficients (Coeff[0], Coeff[1],Coeff[2],Coeff[3]);
933 Coeff[3] = Coeff[3] - Coeff[0]*xloc - Coeff[1]*yloc - Coeff[2]*zloc;
934 return myEpsilon = Compute(S,Standard_False,Coeff,loc,dim,g,inertia,Eps);
935}
936
937Standard_Real GProp_VGProps::Perform(Face& S, Domain& D, const gp_Pln& Pl, const Standard_Real Eps){
938 Standard_Real xloc, yloc, zloc;
939 loc.Coord (xloc, yloc, zloc);
940 Standard_Real Coeff[4];
941 Pl.Coefficients (Coeff[0], Coeff[1],Coeff[2],Coeff[3]);
942 Coeff[3] = Coeff[3] - Coeff[0]*xloc - Coeff[1]*yloc - Coeff[2]*zloc;
943 return myEpsilon = Compute(S,D,Standard_False,Coeff,loc,dim,g,inertia,Eps);
944}
945
946void GProp_VGProps::Perform(const Face& S, const gp_Pln& Pl){
947 Standard_Real xloc, yloc, zloc;
948 loc.Coord (xloc, yloc, zloc);
949 Standard_Real Coeff[4];
950 Pl.Coefficients (Coeff[0], Coeff[1],Coeff[2],Coeff[3]);
951 Coeff[3] = Coeff[3] - Coeff[0]*xloc - Coeff[1]*yloc - Coeff[2]*zloc;
952 Compute(S,Standard_False,Coeff,loc,dim,g,inertia);
953 myEpsilon = 1.0;
954 return;
955}
956
957void GProp_VGProps::Perform(Face& S, Domain& D, const gp_Pln& Pl){
958 Standard_Real xloc, yloc, zloc;
959 loc.Coord (xloc, yloc, zloc);
960 Standard_Real Coeff[4];
961 Pl.Coefficients (Coeff[0], Coeff[1],Coeff[2],Coeff[3]);
962 Coeff[3] = Coeff[3] - Coeff[0]*xloc - Coeff[1]*yloc - Coeff[2]*zloc;
963 Compute(S,D,Standard_False,Coeff,loc,dim,g,inertia);
964 myEpsilon = 1.0;
965 return;
966}
967
968Standard_Real GProp_VGProps::GetEpsilon(){
969 return myEpsilon;
970}