0024002: Overall code and build procedure refactoring -- automatic
[occt.git] / src / BRepGProp / BRepGProp_VinertGK.cxx
CommitLineData
b311480e 1// Created on: 2005-12-09
2// Created by: Sergey KHROMOV
973c2be1 3// Copyright (c) 2005-2014 OPEN CASCADE SAS
b311480e 4//
973c2be1 5// This file is part of Open CASCADE Technology software library.
b311480e 6//
d5f74e42 7// This library is free software; you can redistribute it and/or modify it under
8// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1 9// by the Free Software Foundation, with special exception defined in the file
10// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
11// distribution for complete text of the license and disclaimer of any warranty.
b311480e 12//
973c2be1 13// Alternatively, this file may be used under the terms of Open CASCADE
14// commercial license or contractual agreement.
7fd59977 15
424cd6bb 16
42cf5bc1 17#include <BRepGProp_Domain.hxx>
18#include <BRepGProp_Face.hxx>
424cd6bb 19#include <BRepGProp_TFunction.hxx>
42cf5bc1 20#include <BRepGProp_VinertGK.hxx>
21#include <gp_Pln.hxx>
22#include <gp_Pnt.hxx>
23#include <math_KronrodSingleIntegration.hxx>
24#include <TColStd_Array1OfBoolean.hxx>
25#include <TColStd_HArray1OfReal.hxx>
7fd59977 26
27//==========================================================================
28//function : Constructor
7fd59977 29//==========================================================================
424cd6bb 30BRepGProp_VinertGK::BRepGProp_VinertGK(BRepGProp_Face &theSurface,
31 const gp_Pnt &theLocation,
32 const Standard_Real theTolerance,
33 const Standard_Boolean theCGFlag,
34 const Standard_Boolean theIFlag):
35 myErrorReached(0.)
7fd59977 36{
37 SetLocation(theLocation);
38 Perform(theSurface, theTolerance, theCGFlag, theIFlag);
39}
40
41//==========================================================================
42//function : Constructor
43//
44//==========================================================================
45
424cd6bb 46BRepGProp_VinertGK::BRepGProp_VinertGK( BRepGProp_Face &theSurface,
47 const gp_Pnt &thePoint,
48 const gp_Pnt &theLocation,
49 const Standard_Real theTolerance,
50 const Standard_Boolean theCGFlag,
51 const Standard_Boolean theIFlag):
52 myErrorReached(0.)
7fd59977 53{
54 SetLocation(theLocation);
55 Perform(theSurface, thePoint, theTolerance, theCGFlag, theIFlag);
56}
57
58//==========================================================================
59//function : Constructor
60//
61//==========================================================================
62
424cd6bb 63BRepGProp_VinertGK::BRepGProp_VinertGK(BRepGProp_Face &theSurface,
64 BRepGProp_Domain &theDomain,
65 const gp_Pnt &theLocation,
66 const Standard_Real theTolerance,
67 const Standard_Boolean theCGFlag,
68 const Standard_Boolean theIFlag):
69 myErrorReached(0.)
7fd59977 70{
71 SetLocation(theLocation);
72 Perform(theSurface, theDomain, theTolerance, theCGFlag, theIFlag);
73}
74
75//==========================================================================
76//function : Constructor
77//
78//==========================================================================
79
424cd6bb 80BRepGProp_VinertGK::BRepGProp_VinertGK(BRepGProp_Face &theSurface,
81 BRepGProp_Domain &theDomain,
82 const gp_Pnt &thePoint,
83 const gp_Pnt &theLocation,
84 const Standard_Real theTolerance,
85 const Standard_Boolean theCGFlag,
86 const Standard_Boolean theIFlag):
87 myErrorReached(0.)
7fd59977 88{
89 SetLocation(theLocation);
90 Perform(theSurface, theDomain, thePoint, theTolerance, theCGFlag, theIFlag);
91}
92
93//==========================================================================
94//function : Constructor
95//
96//==========================================================================
97
424cd6bb 98BRepGProp_VinertGK::BRepGProp_VinertGK(BRepGProp_Face &theSurface,
99 const gp_Pln &thePlane,
100 const gp_Pnt &theLocation,
101 const Standard_Real theTolerance,
102 const Standard_Boolean theCGFlag,
103 const Standard_Boolean theIFlag):
104 myErrorReached(0.)
7fd59977 105{
106 SetLocation(theLocation);
107 Perform(theSurface, thePlane, theTolerance, theCGFlag, theIFlag);
108}
109
110//==========================================================================
111//function : Constructor
112//
113//==========================================================================
114
424cd6bb 115BRepGProp_VinertGK::BRepGProp_VinertGK(BRepGProp_Face &theSurface,
116 BRepGProp_Domain &theDomain,
117 const gp_Pln &thePlane,
118 const gp_Pnt &theLocation,
119 const Standard_Real theTolerance,
120 const Standard_Boolean theCGFlag,
121 const Standard_Boolean theIFlag):
122 myErrorReached(0.)
7fd59977 123{
124 SetLocation(theLocation);
125 Perform(theSurface, theDomain, thePlane, theTolerance, theCGFlag, theIFlag);
126}
127
128//==========================================================================
129//function : Perform
130// Compute the properties.
131//==========================================================================
132
424cd6bb 133Standard_Real BRepGProp_VinertGK::Perform(BRepGProp_Face &theSurface,
134 const Standard_Real theTolerance,
135 const Standard_Boolean theCGFlag,
136 const Standard_Boolean theIFlag)
7fd59977 137
138{
139 Standard_Real aShift[] = { 0., 0., 0. };
140
141 return PrivatePerform(theSurface, NULL, Standard_True, &aShift, theTolerance,
424cd6bb 142 theCGFlag, theIFlag);
7fd59977 143}
144
145//==========================================================================
146//function : Perform
147// Compute the properties.
148//==========================================================================
149
424cd6bb 150Standard_Real BRepGProp_VinertGK::Perform(BRepGProp_Face &theSurface,
151 const gp_Pnt &thePoint,
152 const Standard_Real theTolerance,
153 const Standard_Boolean theCGFlag,
154 const Standard_Boolean theIFlag)
7fd59977 155
156{
157 gp_XYZ aXYZ(thePoint.XYZ().Subtracted(loc.XYZ()));
158 Standard_Real aShift[3];
159
160 aXYZ.Coord(aShift[0], aShift[1], aShift[2]);
161
162 return PrivatePerform(theSurface, NULL, Standard_True, &aShift, theTolerance,
424cd6bb 163 theCGFlag, theIFlag);
7fd59977 164}
165
166//==========================================================================
167//function : Perform
168// Compute the properties.
169//==========================================================================
170
424cd6bb 171Standard_Real BRepGProp_VinertGK::Perform(BRepGProp_Face &theSurface,
172 BRepGProp_Domain &theDomain,
173 const Standard_Real theTolerance,
174 const Standard_Boolean theCGFlag,
175 const Standard_Boolean theIFlag)
7fd59977 176
177{
178 Standard_Real aShift[] = { 0., 0., 0. };
179
180 return PrivatePerform(theSurface, &theDomain,
424cd6bb 181 Standard_True, &aShift, theTolerance,
182 theCGFlag, theIFlag);
7fd59977 183}
184
185//==========================================================================
186//function : Perform
187// Compute the properties.
188//==========================================================================
189
424cd6bb 190Standard_Real BRepGProp_VinertGK::Perform(BRepGProp_Face &theSurface,
191 BRepGProp_Domain &theDomain,
192 const gp_Pnt &thePoint,
193 const Standard_Real theTolerance,
194 const Standard_Boolean theCGFlag,
195 const Standard_Boolean theIFlag)
7fd59977 196
197{
198 gp_XYZ aXYZ(thePoint.XYZ().Subtracted(loc.XYZ()));
199 Standard_Real aShift[3];
200
201 aXYZ.Coord(aShift[0], aShift[1], aShift[2]);
202
203 return PrivatePerform(theSurface, &theDomain,
424cd6bb 204 Standard_True, &aShift, theTolerance,
205 theCGFlag, theIFlag);
7fd59977 206}
207
208//==========================================================================
209//function : Perform
210// Compute the properties.
211//==========================================================================
212
424cd6bb 213Standard_Real BRepGProp_VinertGK::Perform(BRepGProp_Face &theSurface,
214 const gp_Pln &thePlane,
215 const Standard_Real theTolerance,
216 const Standard_Boolean theCGFlag,
217 const Standard_Boolean theIFlag)
7fd59977 218
219{
220 Standard_Real aCoeff[4];
221 Standard_Real aXLoc;
222 Standard_Real aYLoc;
223 Standard_Real aZLoc;
224
225 loc.Coord(aXLoc, aYLoc, aZLoc);
226 thePlane.Coefficients (aCoeff[0], aCoeff[1], aCoeff[2], aCoeff[3]);
227 aCoeff[3] = aCoeff[3] - aCoeff[0]*aXLoc - aCoeff[1]*aYLoc - aCoeff[2]*aZLoc;
228
229 return PrivatePerform(theSurface, NULL,
424cd6bb 230 Standard_False, &aCoeff, theTolerance,
231 theCGFlag, theIFlag);
7fd59977 232}
233
234//==========================================================================
235//function : Perform
236// Compute the properties.
237//==========================================================================
238
424cd6bb 239Standard_Real BRepGProp_VinertGK::Perform(BRepGProp_Face &theSurface,
240 BRepGProp_Domain &theDomain,
241 const gp_Pln &thePlane,
242 const Standard_Real theTolerance,
243 const Standard_Boolean theCGFlag,
244 const Standard_Boolean theIFlag)
7fd59977 245
246{
247 Standard_Real aCoeff[4];
248 Standard_Real aXLoc;
249 Standard_Real aYLoc;
250 Standard_Real aZLoc;
251
252 loc.Coord(aXLoc, aYLoc, aZLoc);
253 thePlane.Coefficients (aCoeff[0], aCoeff[1], aCoeff[2], aCoeff[3]);
254 aCoeff[3] = aCoeff[3] - aCoeff[0]*aXLoc - aCoeff[1]*aYLoc - aCoeff[2]*aZLoc;
255
256 return PrivatePerform(theSurface, &theDomain,
424cd6bb 257 Standard_False, &aCoeff, theTolerance,
258 theCGFlag, theIFlag);
7fd59977 259}
260
261//==========================================================================
262//function : PrivatePerform
263// Compute the properties.
264//==========================================================================
265
424cd6bb 266Standard_Real BRepGProp_VinertGK::PrivatePerform
267(BRepGProp_Face &theSurface,
268 const Standard_Address thePtrDomain,
269 const Standard_Boolean IsByPoint,
270 const Standard_Address theCoeffs,
271 const Standard_Real theTolerance,
272 const Standard_Boolean theCGFlag,
273 const Standard_Boolean theIFlag)
7fd59977 274
275{
276
277 const Standard_Real aTTol = 1.e-9;
278 Standard_Real *aCoeffs = (Standard_Real *)theCoeffs;
279
280 // Compute the number of 2d bounding curves of the face.
424cd6bb 281 BRepGProp_Domain *aPDomain = NULL;
7fd59977 282 Standard_Integer aNbCurves = 0;
283
284 // If the pointer to the domain is NULL, there is only one curve to treat:
285 // U isoline with the UMax parameter.
286 if (thePtrDomain == NULL)
287 aNbCurves = 1;
288 else {
424cd6bb 289 aPDomain = (BRepGProp_Domain *)thePtrDomain;
7fd59977 290
291 for (aPDomain->Init(); aPDomain->More(); aPDomain->Next())
292 aNbCurves++;
293 }
294
295 if (aNbCurves == 0) {
296 myErrorReached = -1.;
297
298 return myErrorReached;
299 }
300
301 //Standard_Real aCrvTol = 0.5*theTolerance/aNbCurves;
302 Standard_Real aCrvTol = 0.1*theTolerance;
303 Standard_Real aUMin;
304 Standard_Real aUMax;
305 Standard_Real aTMin;
306 Standard_Real aTMax;
307 Standard_Integer aNbPnts;
308 Standard_Integer aNbMaxIter = 1000;
309 Standard_Integer aNbVal = 10;
310 Standard_Integer k;
311 math_Vector aLocalValue(1, aNbVal);
312 math_Vector aLocalTolReached(1, aNbVal);
313 math_Vector aValue(1, aNbVal);
314 math_Vector aTolReached(1, aNbVal);
315 TColStd_Array1OfBoolean CFlags(1, aNbVal);
316 CFlags.Init(Standard_False);
317 Standard_Boolean isMore;
318
319 //aNbVal = 1;
320 aValue.Init(0.);
321 aTolReached.Init(0.);
322
323 CFlags.Init(Standard_False);
324 CFlags(1) = Standard_True;
325
326 if(theCGFlag || theIFlag) {
327 Standard_Integer i;
328 for(i = 2; i <= 4; ++i) {CFlags(i) = Standard_True;}
329 }
330
331 if(theIFlag) {
332 Standard_Integer i;
333 for(i = 5; i <= 10; ++i) {CFlags(i) = Standard_True;}
334 }
335
336 theSurface.Bounds(aUMin, aUMax, aTMin, aTMax);
337
338 if (thePtrDomain == NULL)
339 isMore = Standard_True;
340 else {
341 aPDomain->Init();
342 isMore = aPDomain->More();
343 }
344
345 while(isMore) {
346 // If the pointer to the domain is NULL, there is only one curve to treat:
347 // U isoline with the UMax parameter.
348
349 if (thePtrDomain == NULL)
350 theSurface.Load(Standard_False, GeomAbs_IsoU);
351 else
352 theSurface.Load(aPDomain->Value());
353
354 aTMin = theSurface.FirstParameter();
355 aTMax = theSurface.LastParameter();
356
357
358 // Get the spans on the curve.
359 Handle(TColStd_HArray1OfReal) aTKnots;
424cd6bb 360 BRepGProp_TFunction aTFunc(theSurface, loc, IsByPoint, theCoeffs,
361 aUMin, aCrvTol);
7fd59977 362
363 theSurface.GetTKnots(aTMin, aTMax, aTKnots);
364
365 Standard_Integer iU = aTKnots->Upper();
366 Standard_Integer aNbTIntervals = aTKnots->Length() - 1;
367 //Standard_Real aTolSpan = aCrvTol/aNbTIntervals;
368 Standard_Real aTolSpan = 0.9*theTolerance; //Relative error
369 math_KronrodSingleIntegration anIntegral;
370 GProp_ValueType aValueType;
371
372
373 // Empirical criterion.
374 aNbPnts = Min(15, theSurface.IntegrationOrder()/aNbTIntervals + 1);
375 aNbPnts = Max(5, aNbPnts);
424cd6bb 376 // aNbPnts = theSurface.IntegrationOrder();
7fd59977 377
378 aLocalValue.Init(0.);
379 aLocalTolReached.Init(0.);
380
381 for (k = 1; k <= aNbVal; k++) {
382
383 if(!CFlags(k)) continue;
384
385 Standard_Integer i = aTKnots->Lower();
386
387 switch (k) {
388 case 1: aValueType = GProp_Mass; break;
389 case 2: aValueType = GProp_CenterMassX; break;
390 case 3: aValueType = GProp_CenterMassY; break;
391 case 4: aValueType = GProp_CenterMassZ; break;
392 case 5: aValueType = GProp_InertiaXX; break;
393 case 6: aValueType = GProp_InertiaYY; break;
394 case 7: aValueType = GProp_InertiaZZ; break;
395 case 8: aValueType = GProp_InertiaXY; break;
396 case 9: aValueType = GProp_InertiaXZ; break;
397 case 10: aValueType = GProp_InertiaYZ; break;
398
399 default: myErrorReached = -1.; return myErrorReached;
400 }
401 aTFunc.SetValueType(aValueType);
402
403 Standard_Real err1 = 0.;
404 while (i < iU) {
405
424cd6bb 406 //cout << "-------------- Span " << i << " nbp: " << aNbPnts << endl;
407 Standard_Real aT1 = aTKnots->Value(i++);
408 Standard_Real aT2 = aTKnots->Value(i);
7fd59977 409
424cd6bb 410 if(aT2 - aT1 < aTTol) continue;
7fd59977 411
424cd6bb 412 aTFunc.SetNbKronrodPoints(aNbPnts);
413 aTFunc.Init();
414 aTFunc.SetTolerance(aCrvTol/(aT2-aT1));
415 anIntegral.Perform(aTFunc, aT1, aT2, aNbPnts, aTolSpan, aNbMaxIter);
7fd59977 416
424cd6bb 417 if (!anIntegral.IsDone()) {
418 myErrorReached = -1.;
7fd59977 419
424cd6bb 420 return myErrorReached;
421 }
7fd59977 422
424cd6bb 423 aLocalValue(k) += anIntegral.Value();
424 err1 = aTFunc.AbsolutError()*(aT2 - aT1);
425 //cout << "Errors: " << anIntegral.NbIterReached() << " " << anIntegral.AbsolutError() << " " << err1 << endl;
426 aLocalTolReached(k) += anIntegral.AbsolutError() + err1;
427 //cout << "--- Errors: " << anIntegral.NbIterReached() << " " << anIntegral.AbsolutError() << " " << err1 << endl;
7fd59977 428 }
429
430 aValue(k) += aLocalValue(k);
431 aTolReached(k) += aLocalTolReached(k);
432 }
433
434 // If the pointer to the domain is NULL, there is only one curve to treat:
435 // U isoline with the UMax parameter.
436 if (thePtrDomain == NULL)
437 isMore = Standard_False;
438 else {
439 aPDomain->Next();
440 isMore = aPDomain->More();
441 }
442 }
443
444 // Get volume value.
445 dim = aValue(1);
446 myErrorReached = aTolReached(1);
447 myAbsolutError = myErrorReached;
448 Standard_Real anAbsDim = Abs(dim);
449 Standard_Real aVolTol = Epsilon(myAbsolutError);
450 if(anAbsDim >= aVolTol) myErrorReached /= anAbsDim;
451
452 if(theCGFlag || theIFlag) {
453 // Compute values of center of mass.
454 if(anAbsDim >= aVolTol) {
455 if (IsByPoint) {
424cd6bb 456 aValue(2) = aCoeffs[0] + aValue(2)/dim;
457 aValue(3) = aCoeffs[1] + aValue(3)/dim;
458 aValue(4) = aCoeffs[2] + aValue(4)/dim;
7fd59977 459 } else {
424cd6bb 460 aValue(2) /= dim;
461 aValue(3) /= dim;
462 aValue(4) /= dim;
7fd59977 463 }
464 } else {
465 aValue(2) = 0.;
466 aValue(3) = 0.;
467 aValue(4) = 0.;
468 dim = 0.;
469 }
470 g.SetCoord(aValue(2), aValue(3), aValue(4));
471 }
472
473 if(theIFlag) {
474 // Fill the matrix of inertia.
475 inertia.SetCols (gp_XYZ (aValue(5), aValue(8), aValue(9)),
424cd6bb 476 gp_XYZ (aValue(8), aValue(6), aValue(10)),
477 gp_XYZ (aValue(9), aValue(10), aValue(7)));
7fd59977 478 }
479 //return myErrorReached;
480 return myAbsolutError;
481}
424cd6bb 482