6b84c3f7 |
1 | // Created on: 2015-09-21 |
2 | // Copyright (c) 2015 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 | #include <GeomEvaluator_OffsetSurface.hxx> |
16 | |
17 | #include <GeomAdaptor_HSurface.hxx> |
18 | #include <CSLib.hxx> |
19 | #include <CSLib_NormalStatus.hxx> |
20 | #include <Geom_BezierSurface.hxx> |
21 | #include <Geom_BSplineSurface.hxx> |
22 | #include <Geom_UndefinedValue.hxx> |
23 | #include <GeomAdaptor_HSurface.hxx> |
24 | #include <gp_Vec.hxx> |
25 | #include <Standard_RangeError.hxx> |
26 | |
27 | |
92efcf78 |
28 | IMPLEMENT_STANDARD_RTTIEXT(GeomEvaluator_OffsetSurface,GeomEvaluator_Surface) |
29 | |
6b84c3f7 |
30 | template<class SurfOrAdapt> |
31 | static void derivatives(Standard_Integer theMaxOrder, |
32 | Standard_Integer theMinOrder, |
33 | const Standard_Real theU, |
34 | const Standard_Real theV, |
35 | const SurfOrAdapt& theBasisSurf, |
36 | const Standard_Integer theNU, |
37 | const Standard_Integer theNV, |
38 | const Standard_Boolean theAlongU, |
39 | const Standard_Boolean theAlongV, |
40 | const Handle(Geom_BSplineSurface)& theL, |
41 | TColgp_Array2OfVec& theDerNUV, |
42 | TColgp_Array2OfVec& theDerSurf); |
43 | |
44 | |
45 | |
46 | GeomEvaluator_OffsetSurface::GeomEvaluator_OffsetSurface( |
47 | const Handle(Geom_Surface)& theBase, |
48 | const Standard_Real theOffset, |
49 | const Handle(Geom_OsculatingSurface)& theOscSurf) |
50 | : GeomEvaluator_Surface(), |
51 | myBaseSurf(theBase), |
52 | myOffset(theOffset), |
53 | myOscSurf(theOscSurf) |
54 | { |
55 | if (!myOscSurf.IsNull()) |
56 | return; // osculating surface already exists |
57 | |
58 | // Create osculating surface for B-spline and Besier surfaces only |
59 | if (myBaseSurf->IsKind(STANDARD_TYPE(Geom_BSplineSurface)) || |
60 | myBaseSurf->IsKind(STANDARD_TYPE(Geom_BezierSurface))) |
61 | myOscSurf = new Geom_OsculatingSurface(myBaseSurf, Precision::Confusion()); |
62 | } |
63 | |
64 | GeomEvaluator_OffsetSurface::GeomEvaluator_OffsetSurface( |
65 | const Handle(GeomAdaptor_HSurface)& theBase, |
66 | const Standard_Real theOffset, |
67 | const Handle(Geom_OsculatingSurface)& theOscSurf) |
68 | : GeomEvaluator_Surface(), |
69 | myBaseAdaptor(theBase), |
525ec87c |
70 | myOffset(theOffset), |
71 | myOscSurf(theOscSurf) |
6b84c3f7 |
72 | { |
6b84c3f7 |
73 | } |
74 | |
75 | void GeomEvaluator_OffsetSurface::D0( |
76 | const Standard_Real theU, const Standard_Real theV, gp_Pnt& theValue) const |
77 | { |
78 | gp_Vec aD1U, aD1V; |
79 | BaseD1(theU, theV, theValue, aD1U, aD1V); |
80 | CalculateD0(theU, theV, theValue, aD1U, aD1V); |
81 | } |
82 | |
83 | void GeomEvaluator_OffsetSurface::D1( |
84 | const Standard_Real theU, const Standard_Real theV, |
85 | gp_Pnt& theValue, gp_Vec& theD1U, gp_Vec& theD1V) const |
86 | { |
87 | gp_Vec aD2U, aD2V, aD2UV; |
88 | BaseD2(theU, theV, theValue, theD1U, theD1V, aD2U, aD2V, aD2UV); |
89 | CalculateD1(theU, theV, theValue, theD1U, theD1V, aD2U, aD2V, aD2UV); |
90 | } |
91 | |
92 | void GeomEvaluator_OffsetSurface::D2( |
93 | const Standard_Real theU, const Standard_Real theV, |
94 | gp_Pnt& theValue, gp_Vec& theD1U, gp_Vec& theD1V, |
95 | gp_Vec& theD2U, gp_Vec& theD2V, gp_Vec& theD2UV) const |
96 | { |
97 | gp_Vec aD3U, aD3V, aD3UUV, aD3UVV; |
98 | BaseD3(theU, theV, theValue, theD1U, theD1V, |
99 | theD2U, theD2V, theD2UV, aD3U, aD3V, aD3UUV, aD3UVV); |
100 | CalculateD2(theU, theV, theValue, theD1U, theD1V, |
101 | theD2U, theD2V, theD2UV, aD3U, aD3V, aD3UUV, aD3UVV); |
102 | } |
103 | |
104 | void GeomEvaluator_OffsetSurface::D3( |
105 | const Standard_Real theU, const Standard_Real theV, |
106 | gp_Pnt& theValue, gp_Vec& theD1U, gp_Vec& theD1V, |
107 | gp_Vec& theD2U, gp_Vec& theD2V, gp_Vec& theD2UV, |
108 | gp_Vec& theD3U, gp_Vec& theD3V, gp_Vec& theD3UUV, gp_Vec& theD3UVV) const |
109 | { |
110 | BaseD3(theU, theV, theValue, theD1U, theD1V, |
111 | theD2U, theD2V, theD2UV, theD3U, theD3V, theD3UUV, theD3UVV); |
112 | CalculateD3(theU, theV, theValue, theD1U, theD1V, |
113 | theD2U, theD2V, theD2UV, theD3U, theD3V, theD3UUV, theD3UVV); |
114 | } |
115 | |
116 | gp_Vec GeomEvaluator_OffsetSurface::DN( |
117 | const Standard_Real theU, const Standard_Real theV, |
118 | const Standard_Integer theDerU, const Standard_Integer theDerV) const |
119 | { |
120 | Standard_RangeError_Raise_if(theDerU < 0, "GeomEvaluator_OffsetSurface::DN(): theDerU < 0"); |
121 | Standard_RangeError_Raise_if(theDerV < 0, "GeomEvaluator_OffsetSurface::DN(): theDerV < 0"); |
122 | Standard_RangeError_Raise_if(theDerU + theDerV < 1, |
123 | "GeomEvaluator_OffsetSurface::DN(): theDerU + theDerV < 1"); |
124 | |
125 | gp_Pnt aP; |
126 | gp_Vec aD1U, aD1V; |
127 | BaseD1(theU, theV, aP, aD1U, aD1V); |
128 | return CalculateDN(theU, theV, theDerU, theDerV, aD1U, aD1V); |
129 | } |
130 | |
131 | |
132 | void GeomEvaluator_OffsetSurface::BaseD0(const Standard_Real theU, const Standard_Real theV, |
133 | gp_Pnt& theValue) const |
134 | { |
135 | if (!myBaseAdaptor.IsNull()) |
136 | myBaseAdaptor->D0(theU, theV, theValue); |
137 | else |
138 | myBaseSurf->D0(theU, theV, theValue); |
139 | } |
140 | |
141 | void GeomEvaluator_OffsetSurface::BaseD1(const Standard_Real theU, const Standard_Real theV, |
142 | gp_Pnt& theValue, gp_Vec& theD1U, gp_Vec& theD1V) const |
143 | { |
144 | if (!myBaseAdaptor.IsNull()) |
145 | myBaseAdaptor->D1(theU, theV, theValue, theD1U, theD1V); |
146 | else |
147 | myBaseSurf->D1(theU, theV, theValue, theD1U, theD1V); |
148 | } |
149 | |
150 | void GeomEvaluator_OffsetSurface::BaseD2(const Standard_Real theU, const Standard_Real theV, |
151 | gp_Pnt& theValue, gp_Vec& theD1U, gp_Vec& theD1V, |
152 | gp_Vec& theD2U, gp_Vec& theD2V, gp_Vec& theD2UV) const |
153 | { |
154 | if (!myBaseAdaptor.IsNull()) |
155 | myBaseAdaptor->D2(theU, theV, theValue, theD1U, theD1V, theD2U, theD2V, theD2UV); |
156 | else |
157 | myBaseSurf->D2(theU, theV, theValue, theD1U, theD1V, theD2U, theD2V, theD2UV); |
158 | } |
159 | |
160 | void GeomEvaluator_OffsetSurface::BaseD3( |
161 | const Standard_Real theU, const Standard_Real theV, |
162 | gp_Pnt& theValue, gp_Vec& theD1U, gp_Vec& theD1V, |
163 | gp_Vec& theD2U, gp_Vec& theD2V, gp_Vec& theD2UV, |
164 | gp_Vec& theD3U, gp_Vec& theD3V, gp_Vec& theD3UUV, gp_Vec& theD3UVV) const |
165 | { |
166 | if (!myBaseAdaptor.IsNull()) |
167 | myBaseAdaptor->D3(theU, theV, theValue, theD1U, theD1V, |
168 | theD2U, theD2V, theD2UV, theD3U, theD3V, theD3UUV, theD3UVV); |
169 | else |
170 | myBaseSurf->D3(theU, theV, theValue, theD1U, theD1V, |
171 | theD2U, theD2V, theD2UV, theD3U, theD3V, theD3UUV, theD3UVV); |
172 | } |
173 | |
174 | |
175 | void GeomEvaluator_OffsetSurface::CalculateD0( |
176 | const Standard_Real theU, const Standard_Real theV, |
177 | gp_Pnt& theValue, |
178 | const gp_Vec& theD1U, const gp_Vec& theD1V) const |
179 | { |
180 | const Standard_Real MagTol = 1.e-9; |
181 | |
182 | // Normalize derivatives before normal calculation because it gives more stable result. |
183 | // There will be normalized only derivatives greater than 1.0 to avoid differences in last significant digit |
184 | gp_Vec aD1U(theD1U); |
185 | gp_Vec aD1V(theD1V); |
186 | Standard_Real aD1UNorm2 = aD1U.SquareMagnitude(); |
187 | Standard_Real aD1VNorm2 = aD1V.SquareMagnitude(); |
188 | if (aD1UNorm2 > 1.0) |
189 | aD1U /= Sqrt(aD1UNorm2); |
190 | if (aD1VNorm2 > 1.0) |
191 | aD1V /= Sqrt(aD1VNorm2); |
192 | |
193 | gp_Vec aNorm = aD1U.Crossed(aD1V); |
194 | if (aNorm.SquareMagnitude() > MagTol * MagTol) |
195 | { |
196 | // Non singular case. Simple computations. |
197 | aNorm.Normalize(); |
198 | theValue.SetXYZ(theValue.XYZ() + myOffset * aNorm.XYZ()); |
199 | } |
200 | else |
201 | { |
202 | const Standard_Integer MaxOrder = 3; |
203 | |
204 | Handle(Geom_BSplineSurface) L; |
205 | Standard_Boolean isOpposite = Standard_False; |
206 | Standard_Boolean AlongU = Standard_False; |
207 | Standard_Boolean AlongV = Standard_False; |
208 | if (!myOscSurf.IsNull()) |
209 | { |
210 | AlongU = myOscSurf->UOscSurf(theU, theV, isOpposite, L); |
211 | AlongV = myOscSurf->VOscSurf(theU, theV, isOpposite, L); |
212 | } |
213 | const Standard_Real aSign = ((AlongV || AlongU) && isOpposite) ? -1. : 1.; |
214 | |
215 | TColgp_Array2OfVec DerNUV(0, MaxOrder, 0, MaxOrder); |
216 | TColgp_Array2OfVec DerSurf(0, MaxOrder + 1, 0, MaxOrder + 1); |
217 | Standard_Integer OrderU, OrderV; |
218 | Standard_Real Umin = 0, Umax = 0, Vmin = 0, Vmax = 0; |
219 | Bounds(Umin, Umax, Vmin, Vmax); |
220 | |
221 | DerSurf.SetValue(1, 0, theD1U); |
222 | DerSurf.SetValue(0, 1, theD1V); |
223 | if (!myBaseSurf.IsNull()) |
224 | derivatives(MaxOrder, 1, theU, theV, myBaseSurf, 0, 0, AlongU, AlongV, L, DerNUV, DerSurf); |
225 | else |
226 | derivatives(MaxOrder, 1, theU, theV, myBaseAdaptor, 0, 0, AlongU, AlongV, L, DerNUV, DerSurf); |
227 | |
228 | gp_Dir Normal; |
229 | CSLib_NormalStatus NStatus; |
230 | CSLib::Normal(MaxOrder, DerNUV, MagTol, theU, theV, Umin, Umax, Vmin, Vmax, |
231 | NStatus, Normal, OrderU, OrderV); |
232 | if (NStatus == CSLib_InfinityOfSolutions) |
233 | { |
234 | // Replace zero derivative and try to calculate normal |
235 | gp_Vec aNewDU = theD1U; |
236 | gp_Vec aNewDV = theD1V; |
237 | if (ReplaceDerivative(theU, theV, aNewDU, aNewDV, MagTol * MagTol)) |
238 | CSLib::Normal(aNewDU, aNewDV, MagTol, NStatus, Normal); |
239 | } |
240 | |
241 | if (NStatus != CSLib_Defined) |
242 | Geom_UndefinedValue::Raise( |
243 | "GeomEvaluator_OffsetSurface::CalculateD0(): Unable to calculate normal"); |
244 | |
245 | theValue.SetXYZ(theValue.XYZ() + myOffset * aSign * Normal.XYZ()); |
246 | } |
247 | } |
248 | |
249 | void GeomEvaluator_OffsetSurface::CalculateD1( |
250 | const Standard_Real theU, const Standard_Real theV, |
251 | gp_Pnt& theValue, gp_Vec& theD1U, gp_Vec& theD1V, |
252 | const gp_Vec& theD2U, const gp_Vec& theD2V, const gp_Vec& theD2UV) const |
253 | { |
254 | const Standard_Real MagTol = 1.e-9; |
255 | |
256 | // Check offset side. |
257 | Handle(Geom_BSplineSurface) L; |
258 | Standard_Boolean isOpposite = Standard_False; |
259 | Standard_Boolean AlongU = Standard_False; |
260 | Standard_Boolean AlongV = Standard_False; |
6b84c3f7 |
261 | |
262 | // Normalize derivatives before normal calculation because it gives more stable result. |
263 | // There will be normalized only derivatives greater than 1.0 to avoid differences in last significant digit |
264 | gp_Vec aD1U(theD1U); |
265 | gp_Vec aD1V(theD1V); |
266 | Standard_Real aD1UNorm2 = aD1U.SquareMagnitude(); |
267 | Standard_Real aD1VNorm2 = aD1V.SquareMagnitude(); |
268 | if (aD1UNorm2 > 1.0) |
269 | aD1U /= Sqrt(aD1UNorm2); |
270 | if (aD1VNorm2 > 1.0) |
271 | aD1V /= Sqrt(aD1VNorm2); |
272 | |
95ae6ebb |
273 | Standard_Boolean isSingular = Standard_False; |
6b84c3f7 |
274 | Standard_Integer MaxOrder = 0; |
275 | gp_Vec aNorm = aD1U.Crossed(aD1V); |
95ae6ebb |
276 | if (aNorm.SquareMagnitude() <= MagTol * MagTol) |
6b84c3f7 |
277 | { |
95ae6ebb |
278 | if (!myOscSurf.IsNull()) |
6b84c3f7 |
279 | { |
95ae6ebb |
280 | AlongU = myOscSurf->UOscSurf(theU, theV, isOpposite, L); |
281 | AlongV = myOscSurf->VOscSurf(theU, theV, isOpposite, L); |
6b84c3f7 |
282 | } |
95ae6ebb |
283 | |
6b84c3f7 |
284 | MaxOrder = 3; |
95ae6ebb |
285 | isSingular = Standard_True; |
286 | } |
6b84c3f7 |
287 | |
95ae6ebb |
288 | const Standard_Real aSign = ((AlongV || AlongU) && isOpposite) ? -1. : 1.; |
289 | |
290 | if (!isSingular) |
291 | { |
292 | // AlongU or AlongV leads to more complex D1 computation |
293 | // Try to compute D0 and D1 much simpler |
294 | aNorm.Normalize(); |
295 | theValue.SetXYZ(theValue.XYZ() + myOffset * aSign * aNorm.XYZ()); |
296 | |
297 | gp_Vec aN0(aNorm.XYZ()), aN1U, aN1V; |
298 | Standard_Real aScale = (theD1U^theD1V).Dot(aN0); |
299 | aN1U.SetX(theD2U.Y() * theD1V.Z() + theD1U.Y() * theD2UV.Z() |
300 | - theD2U.Z() * theD1V.Y() - theD1U.Z() * theD2UV.Y()); |
301 | aN1U.SetY((theD2U.X() * theD1V.Z() + theD1U.X() * theD2UV.Z() |
302 | - theD2U.Z() * theD1V.X() - theD1U.Z() * theD2UV.X()) * -1.0); |
303 | aN1U.SetZ(theD2U.X() * theD1V.Y() + theD1U.X() * theD2UV.Y() |
304 | - theD2U.Y() * theD1V.X() - theD1U.Y() * theD2UV.X()); |
305 | Standard_Real aScaleU = aN1U.Dot(aN0); |
306 | aN1U.Subtract(aScaleU * aN0); |
307 | aN1U /= aScale; |
308 | |
309 | aN1V.SetX(theD2UV.Y() * theD1V.Z() + theD2V.Z() * theD1U.Y() |
310 | - theD2UV.Z() * theD1V.Y() - theD2V.Y() * theD1U.Z()); |
311 | aN1V.SetY((theD2UV.X() * theD1V.Z() + theD2V.Z() * theD1U.X() |
312 | - theD2UV.Z() * theD1V.X() - theD2V.X() * theD1U.Z()) * -1.0); |
313 | aN1V.SetZ(theD2UV.X() * theD1V.Y() + theD2V.Y() * theD1U.X() |
314 | - theD2UV.Y() * theD1V.X() - theD2V.X() * theD1U.Y()); |
315 | Standard_Real aScaleV = aN1V.Dot(aN0); |
316 | aN1V.Subtract(aScaleV * aN0); |
317 | aN1V /= aScale; |
318 | |
319 | theD1U += myOffset * aSign * aN1U; |
320 | theD1V += myOffset * aSign * aN1V; |
321 | |
322 | return; |
323 | } |
324 | |
6b84c3f7 |
325 | Standard_Integer OrderU, OrderV; |
326 | TColgp_Array2OfVec DerNUV(0, MaxOrder + 1, 0, MaxOrder + 1); |
327 | TColgp_Array2OfVec DerSurf(0, MaxOrder + 2, 0, MaxOrder + 2); |
328 | Standard_Real Umin = 0, Umax = 0, Vmin = 0, Vmax = 0; |
329 | Bounds(Umin, Umax, Vmin, Vmax); |
330 | |
331 | DerSurf.SetValue(1, 0, theD1U); |
332 | DerSurf.SetValue(0, 1, theD1V); |
333 | DerSurf.SetValue(1, 1, theD2UV); |
334 | DerSurf.SetValue(2, 0, theD2U); |
335 | DerSurf.SetValue(0, 2, theD2V); |
336 | if (!myBaseSurf.IsNull()) |
337 | derivatives(MaxOrder, 2, theU, theV, myBaseSurf, 1, 1, AlongU, AlongV, L, DerNUV, DerSurf); |
338 | else |
339 | derivatives(MaxOrder, 2, theU, theV, myBaseAdaptor, 1, 1, AlongU, AlongV, L, DerNUV, DerSurf); |
340 | |
341 | gp_Dir Normal; |
342 | CSLib_NormalStatus NStatus; |
343 | CSLib::Normal(MaxOrder, DerNUV, MagTol, theU, theV, Umin, Umax, Vmin, Vmax, NStatus, Normal, OrderU, OrderV); |
344 | if (NStatus == CSLib_InfinityOfSolutions) |
345 | { |
346 | gp_Vec aNewDU = theD1U; |
347 | gp_Vec aNewDV = theD1V; |
348 | // Replace zero derivative and try to calculate normal |
349 | if (ReplaceDerivative(theU, theV, aNewDU, aNewDV, MagTol * MagTol)) |
350 | { |
351 | DerSurf.SetValue(1, 0, aNewDU); |
352 | DerSurf.SetValue(0, 1, aNewDV); |
353 | if (!myBaseSurf.IsNull()) |
354 | derivatives(MaxOrder, 2, theU, theV, myBaseSurf, 1, 1, AlongU, AlongV, L, DerNUV, DerSurf); |
355 | else |
356 | derivatives(MaxOrder, 2, theU, theV, myBaseAdaptor, 1, 1, AlongU, AlongV, L, DerNUV, DerSurf); |
357 | CSLib::Normal(MaxOrder, DerNUV, MagTol, theU, theV, Umin, Umax, Vmin, Vmax, |
358 | NStatus, Normal, OrderU, OrderV); |
359 | } |
360 | } |
361 | |
362 | if (NStatus != CSLib_Defined) |
363 | Geom_UndefinedValue::Raise( |
364 | "GeomEvaluator_OffsetSurface::CalculateD1(): Unable to calculate normal"); |
365 | |
366 | theValue.SetXYZ(theValue.XYZ() + myOffset * aSign * Normal.XYZ()); |
367 | |
368 | theD1U = DerSurf(1, 0) + myOffset * aSign * CSLib::DNNormal(1, 0, DerNUV, OrderU, OrderV); |
369 | theD1V = DerSurf(0, 1) + myOffset * aSign * CSLib::DNNormal(0, 1, DerNUV, OrderU, OrderV); |
370 | } |
371 | |
372 | void GeomEvaluator_OffsetSurface::CalculateD2( |
373 | const Standard_Real theU, const Standard_Real theV, |
374 | gp_Pnt& theValue, gp_Vec& theD1U, gp_Vec& theD1V, |
375 | gp_Vec& theD2U, gp_Vec& theD2V, gp_Vec& theD2UV, |
376 | const gp_Vec& theD3U, const gp_Vec& theD3V, const gp_Vec& theD3UUV, const gp_Vec& theD3UVV) const |
377 | { |
378 | const Standard_Real MagTol = 1.e-9; |
379 | |
380 | gp_Dir Normal; |
381 | CSLib_NormalStatus NStatus; |
382 | CSLib::Normal(theD1U, theD1V, MagTol, NStatus, Normal); |
383 | |
384 | const Standard_Integer MaxOrder = (NStatus == CSLib_Defined) ? 0 : 3; |
385 | Standard_Integer OrderU, OrderV; |
386 | TColgp_Array2OfVec DerNUV(0, MaxOrder + 2, 0, MaxOrder + 2); |
387 | TColgp_Array2OfVec DerSurf(0, MaxOrder + 3, 0, MaxOrder + 3); |
388 | |
389 | Standard_Real Umin = 0, Umax = 0, Vmin = 0, Vmax = 0; |
390 | Bounds(Umin, Umax, Vmin, Vmax); |
391 | |
392 | DerSurf.SetValue(1, 0, theD1U); |
393 | DerSurf.SetValue(0, 1, theD1V); |
394 | DerSurf.SetValue(1, 1, theD2UV); |
395 | DerSurf.SetValue(2, 0, theD2U); |
396 | DerSurf.SetValue(0, 2, theD2V); |
397 | DerSurf.SetValue(3, 0, theD3U); |
398 | DerSurf.SetValue(2, 1, theD3UUV); |
399 | DerSurf.SetValue(1, 2, theD3UVV); |
400 | DerSurf.SetValue(0, 3, theD3V); |
401 | //********************* |
402 | |
403 | Handle(Geom_BSplineSurface) L; |
404 | Standard_Boolean isOpposite = Standard_False; |
405 | Standard_Boolean AlongU = Standard_False; |
406 | Standard_Boolean AlongV = Standard_False; |
95ae6ebb |
407 | if ((NStatus != CSLib_Defined) && !myOscSurf.IsNull()) |
6b84c3f7 |
408 | { |
409 | AlongU = myOscSurf->UOscSurf(theU, theV, isOpposite, L); |
410 | AlongV = myOscSurf->VOscSurf(theU, theV, isOpposite, L); |
411 | } |
412 | const Standard_Real aSign = ((AlongV || AlongU) && isOpposite) ? -1. : 1.; |
413 | |
414 | if (!myBaseSurf.IsNull()) |
415 | derivatives(MaxOrder, 3, theU, theV, myBaseSurf, 2, 2, AlongU, AlongV, L, DerNUV, DerSurf); |
416 | else |
417 | derivatives(MaxOrder, 3, theU, theV, myBaseAdaptor, 2, 2, AlongU, AlongV, L, DerNUV, DerSurf); |
418 | |
419 | CSLib::Normal(MaxOrder, DerNUV, MagTol, theU, theV, Umin, Umax, Vmin, Vmax, |
420 | NStatus, Normal, OrderU, OrderV); |
421 | if (NStatus != CSLib_Defined) |
422 | Geom_UndefinedValue::Raise( |
423 | "GeomEvaluator_OffsetSurface::CalculateD2(): Unable to calculate normal"); |
424 | |
425 | theValue.SetXYZ(theValue.XYZ() + myOffset * aSign * Normal.XYZ()); |
426 | |
427 | theD1U = DerSurf(1, 0) + myOffset * aSign * CSLib::DNNormal(1, 0, DerNUV, OrderU, OrderV); |
428 | theD1V = DerSurf(0, 1) + myOffset * aSign * CSLib::DNNormal(0, 1, DerNUV, OrderU, OrderV); |
429 | |
430 | if (!myBaseSurf.IsNull()) |
431 | { |
432 | theD2U = myBaseSurf->DN(theU, theV, 2, 0); |
433 | theD2V = myBaseSurf->DN(theU, theV, 0, 2); |
434 | theD2UV = myBaseSurf->DN(theU, theV, 1, 1); |
435 | } |
436 | else |
437 | { |
438 | theD2U = myBaseAdaptor->DN(theU, theV, 2, 0); |
439 | theD2V = myBaseAdaptor->DN(theU, theV, 0, 2); |
440 | theD2UV = myBaseAdaptor->DN(theU, theV, 1, 1); |
441 | } |
442 | |
443 | theD2U += aSign * myOffset * CSLib::DNNormal(2, 0, DerNUV, OrderU, OrderV); |
444 | theD2V += aSign * myOffset * CSLib::DNNormal(0, 2, DerNUV, OrderU, OrderV); |
445 | theD2UV += aSign * myOffset * CSLib::DNNormal(1, 1, DerNUV, OrderU, OrderV); |
446 | } |
447 | |
448 | void GeomEvaluator_OffsetSurface::CalculateD3( |
449 | const Standard_Real theU, const Standard_Real theV, |
450 | gp_Pnt& theValue, gp_Vec& theD1U, gp_Vec& theD1V, |
451 | gp_Vec& theD2U, gp_Vec& theD2V, gp_Vec& theD2UV, |
452 | gp_Vec& theD3U, gp_Vec& theD3V, gp_Vec& theD3UUV, gp_Vec& theD3UVV) const |
453 | { |
454 | const Standard_Real MagTol = 1.e-9; |
455 | |
456 | gp_Dir Normal; |
457 | CSLib_NormalStatus NStatus; |
458 | CSLib::Normal(theD1U, theD1V, MagTol, NStatus, Normal); |
459 | const Standard_Integer MaxOrder = (NStatus == CSLib_Defined) ? 0 : 3; |
460 | Standard_Integer OrderU, OrderV; |
461 | TColgp_Array2OfVec DerNUV(0, MaxOrder + 3, 0, MaxOrder + 3); |
462 | TColgp_Array2OfVec DerSurf(0, MaxOrder + 4, 0, MaxOrder + 4); |
463 | Standard_Real Umin = 0, Umax = 0, Vmin = 0, Vmax = 0; |
464 | Bounds(Umin, Umax, Vmin, Vmax); |
465 | |
466 | DerSurf.SetValue(1, 0, theD1U); |
467 | DerSurf.SetValue(0, 1, theD1V); |
468 | DerSurf.SetValue(1, 1, theD2UV); |
469 | DerSurf.SetValue(2, 0, theD2U); |
470 | DerSurf.SetValue(0, 2, theD2V); |
471 | DerSurf.SetValue(3, 0, theD3U); |
472 | DerSurf.SetValue(2, 1, theD3UUV); |
473 | DerSurf.SetValue(1, 2, theD3UVV); |
474 | DerSurf.SetValue(0, 3, theD3V); |
475 | |
476 | |
477 | //********************* |
478 | Handle(Geom_BSplineSurface) L; |
479 | Standard_Boolean isOpposite = Standard_False; |
480 | Standard_Boolean AlongU = Standard_False; |
481 | Standard_Boolean AlongV = Standard_False; |
95ae6ebb |
482 | if ((NStatus != CSLib_Defined) && !myOscSurf.IsNull()) |
6b84c3f7 |
483 | { |
484 | AlongU = myOscSurf->UOscSurf(theU, theV, isOpposite, L); |
485 | AlongV = myOscSurf->VOscSurf(theU, theV, isOpposite, L); |
486 | } |
487 | const Standard_Real aSign = ((AlongV || AlongU) && isOpposite) ? -1. : 1.; |
488 | |
489 | if (!myBaseSurf.IsNull()) |
490 | derivatives(MaxOrder, 3, theU, theV, myBaseSurf, 3, 3, AlongU, AlongV, L, DerNUV, DerSurf); |
491 | else |
492 | derivatives(MaxOrder, 3, theU, theV, myBaseAdaptor, 3, 3, AlongU, AlongV, L, DerNUV, DerSurf); |
493 | |
494 | CSLib::Normal(MaxOrder, DerNUV, MagTol, theU, theV, Umin, Umax, Vmin, Vmax, |
495 | NStatus, Normal, OrderU, OrderV); |
496 | if (NStatus != CSLib_Defined) |
497 | Geom_UndefinedValue::Raise( |
498 | "GeomEvaluator_OffsetSurface::CalculateD3(): Unable to calculate normal"); |
499 | |
500 | theValue.SetXYZ(theValue.XYZ() + myOffset * aSign * Normal.XYZ()); |
501 | |
502 | theD1U = DerSurf(1, 0) + myOffset * aSign * CSLib::DNNormal(1, 0, DerNUV, OrderU, OrderV); |
503 | theD1V = DerSurf(0, 1) + myOffset * aSign * CSLib::DNNormal(0, 1, DerNUV, OrderU, OrderV); |
504 | |
505 | if (!myBaseSurf.IsNull()) |
506 | { |
507 | theD2U = myBaseSurf->DN(theU, theV, 2, 0); |
508 | theD2V = myBaseSurf->DN(theU, theV, 0, 2); |
509 | theD2UV = myBaseSurf->DN(theU, theV, 1, 1); |
510 | theD3U = myBaseSurf->DN(theU, theV, 3, 0); |
511 | theD3V = myBaseSurf->DN(theU, theV, 0, 3); |
512 | theD3UUV = myBaseSurf->DN(theU, theV, 2, 1); |
513 | theD3UVV = myBaseSurf->DN(theU, theV, 1, 2); |
514 | } |
515 | else |
516 | { |
517 | theD2U = myBaseAdaptor->DN(theU, theV, 2, 0); |
518 | theD2V = myBaseAdaptor->DN(theU, theV, 0, 2); |
519 | theD2UV = myBaseAdaptor->DN(theU, theV, 1, 1); |
520 | theD3U = myBaseAdaptor->DN(theU, theV, 3, 0); |
521 | theD3V = myBaseAdaptor->DN(theU, theV, 0, 3); |
522 | theD3UUV = myBaseAdaptor->DN(theU, theV, 2, 1); |
523 | theD3UVV = myBaseAdaptor->DN(theU, theV, 1, 2); |
524 | } |
525 | |
526 | theD2U += aSign * myOffset * CSLib::DNNormal(2, 0, DerNUV, OrderU, OrderV); |
527 | theD2V += aSign * myOffset * CSLib::DNNormal(0, 2, DerNUV, OrderU, OrderV); |
528 | theD2UV += aSign * myOffset * CSLib::DNNormal(1, 1, DerNUV, OrderU, OrderV); |
529 | theD3U += aSign * myOffset * CSLib::DNNormal(3, 0, DerNUV, OrderU, OrderV); |
530 | theD3V += aSign * myOffset * CSLib::DNNormal(0, 3, DerNUV, OrderU, OrderV); |
531 | theD3UUV += aSign * myOffset * CSLib::DNNormal(2, 1, DerNUV, OrderU, OrderV); |
532 | theD3UVV += aSign * myOffset * CSLib::DNNormal(1, 2, DerNUV, OrderU, OrderV); |
533 | } |
534 | |
535 | gp_Vec GeomEvaluator_OffsetSurface::CalculateDN( |
536 | const Standard_Real theU, const Standard_Real theV, |
537 | const Standard_Integer theNu, const Standard_Integer theNv, |
538 | const gp_Vec& theD1U, const gp_Vec& theD1V) const |
539 | { |
540 | const Standard_Real MagTol = 1.e-9; |
541 | |
542 | gp_Dir Normal; |
543 | CSLib_NormalStatus NStatus; |
544 | CSLib::Normal(theD1U, theD1V, MagTol, NStatus, Normal); |
545 | const Standard_Integer MaxOrder = (NStatus == CSLib_Defined) ? 0 : 3; |
546 | Standard_Integer OrderU, OrderV; |
547 | TColgp_Array2OfVec DerNUV(0, MaxOrder + theNu, 0, MaxOrder + theNu); |
548 | TColgp_Array2OfVec DerSurf(0, MaxOrder + theNu + 1, 0, MaxOrder + theNv + 1); |
549 | |
550 | Standard_Real Umin = 0, Umax = 0, Vmin = 0, Vmax = 0; |
551 | Bounds(Umin, Umax, Vmin, Vmax); |
552 | |
553 | DerSurf.SetValue(1, 0, theD1U); |
554 | DerSurf.SetValue(0, 1, theD1V); |
555 | |
556 | //********************* |
557 | Handle(Geom_BSplineSurface) L; |
558 | // Is there any osculatingsurface along U or V; |
559 | Standard_Boolean isOpposite = Standard_False; |
560 | Standard_Boolean AlongU = Standard_False; |
561 | Standard_Boolean AlongV = Standard_False; |
95ae6ebb |
562 | if ((NStatus != CSLib_Defined) && !myOscSurf.IsNull()) |
6b84c3f7 |
563 | { |
564 | AlongU = myOscSurf->UOscSurf(theU, theV, isOpposite, L); |
565 | AlongV = myOscSurf->VOscSurf(theU, theV, isOpposite, L); |
566 | } |
567 | const Standard_Real aSign = ((AlongV || AlongU) && isOpposite) ? -1. : 1.; |
568 | |
569 | if (!myBaseSurf.IsNull()) |
570 | derivatives(MaxOrder, 1, theU, theV, myBaseSurf, theNu, theNv, AlongU, AlongV, L, DerNUV, DerSurf); |
571 | else |
572 | derivatives(MaxOrder, 1, theU, theV, myBaseAdaptor, theNu, theNv, AlongU, AlongV, L, DerNUV, DerSurf); |
573 | |
574 | CSLib::Normal(MaxOrder, DerNUV, MagTol, theU, theV, Umin, Umax, Vmin, Vmax, |
575 | NStatus, Normal, OrderU, OrderV); |
576 | if (NStatus != CSLib_Defined) |
577 | Geom_UndefinedValue::Raise( |
578 | "GeomEvaluator_OffsetSurface::CalculateDN(): Unable to calculate normal"); |
579 | |
580 | gp_Vec D; |
581 | if (!myBaseSurf.IsNull()) |
582 | D = myBaseSurf->DN(theU, theV, theNu, theNv); |
583 | else |
584 | D = myBaseAdaptor->DN(theU, theV, theNu, theNv); |
585 | |
586 | D += aSign * myOffset * CSLib::DNNormal(theNu, theNv, DerNUV, OrderU, OrderV); |
587 | return D; |
588 | } |
589 | |
590 | |
591 | void GeomEvaluator_OffsetSurface::Bounds(Standard_Real& theUMin, Standard_Real& theUMax, |
592 | Standard_Real& theVMin, Standard_Real& theVMax) const |
593 | { |
594 | if (!myBaseSurf.IsNull()) |
595 | myBaseSurf->Bounds(theUMin, theUMax, theVMin, theVMax); |
596 | else |
597 | { |
598 | theUMin = myBaseAdaptor->FirstUParameter(); |
599 | theUMax = myBaseAdaptor->LastUParameter(); |
600 | theVMin = myBaseAdaptor->FirstVParameter(); |
601 | theVMax = myBaseAdaptor->LastVParameter(); |
602 | } |
603 | } |
604 | |
605 | |
606 | Standard_Boolean GeomEvaluator_OffsetSurface::ReplaceDerivative( |
607 | const Standard_Real theU, const Standard_Real theV, |
608 | gp_Vec& theDU, gp_Vec& theDV, |
609 | const Standard_Real theSquareTol) const |
610 | { |
611 | Standard_Boolean isReplaceDU = theDU.SquareMagnitude() < theSquareTol; |
612 | Standard_Boolean isReplaceDV = theDV.SquareMagnitude() < theSquareTol; |
613 | Standard_Boolean isReplaced = Standard_False; |
614 | if (isReplaceDU ^ isReplaceDV) |
615 | { |
616 | Standard_Real aU = theU; |
617 | Standard_Real aV = theV; |
618 | Standard_Real aUMin = 0, aUMax = 0, aVMin = 0, aVMax = 0; |
619 | Bounds(aUMin, aUMax, aVMin, aVMax); |
620 | |
621 | // Calculate step along non-zero derivative |
622 | Standard_Real aStep; |
623 | Handle(Adaptor3d_HSurface) aSurfAdapt; |
624 | if (!myBaseAdaptor.IsNull()) |
625 | aSurfAdapt = myBaseAdaptor; |
626 | else |
627 | aSurfAdapt = new GeomAdaptor_HSurface(myBaseSurf); |
628 | if (isReplaceDV) |
629 | { |
630 | aStep = Precision::Confusion() * theDU.Magnitude(); |
631 | if (aStep > aUMax - aUMin) |
632 | aStep = (aUMax - aUMin) / 100.; |
633 | } |
634 | else |
635 | { |
636 | aStep = Precision::Confusion() * theDV.Magnitude(); |
637 | if (aStep > aVMax - aVMin) |
638 | aStep = (aVMax - aVMin) / 100.; |
639 | } |
640 | |
641 | gp_Pnt aP; |
642 | gp_Vec aDU, aDV; |
643 | // Step away from currect parametric coordinates and calculate derivatives once again. |
644 | // Replace zero derivative by the obtained. |
645 | for (Standard_Real aStepSign = -1.0; aStepSign <= 1.0 && !isReplaced; aStepSign += 2.0) |
646 | { |
647 | if (isReplaceDV) |
648 | { |
649 | aU = theU + aStepSign * aStep; |
650 | if (aU < aUMin || aU > aUMax) |
651 | continue; |
652 | } |
653 | else |
654 | { |
655 | aV = theV + aStepSign * aStep; |
656 | if (aV < aVMin || aV > aVMax) |
657 | continue; |
658 | } |
659 | |
660 | BaseD1(aU, aV, aP, aDU, aDV); |
661 | |
662 | if (isReplaceDU && aDU.SquareMagnitude() > theSquareTol) |
663 | { |
664 | theDU = aDU; |
665 | isReplaced = Standard_True; |
666 | } |
667 | if (isReplaceDV && aDV.SquareMagnitude() > theSquareTol) |
668 | { |
669 | theDV = aDV; |
670 | isReplaced = Standard_True; |
671 | } |
672 | } |
673 | } |
674 | return isReplaced; |
675 | } |
676 | |
677 | |
678 | |
679 | |
680 | |
681 | // ===================== Auxiliary functions =================================== |
682 | template<class SurfOrAdapt> |
683 | void derivatives(Standard_Integer theMaxOrder, |
684 | Standard_Integer theMinOrder, |
685 | const Standard_Real theU, |
686 | const Standard_Real theV, |
687 | const SurfOrAdapt& theBasisSurf, |
688 | const Standard_Integer theNU, |
689 | const Standard_Integer theNV, |
690 | const Standard_Boolean theAlongU, |
691 | const Standard_Boolean theAlongV, |
692 | const Handle(Geom_BSplineSurface)& theL, |
693 | TColgp_Array2OfVec& theDerNUV, |
694 | TColgp_Array2OfVec& theDerSurf) |
695 | { |
696 | Standard_Integer i, j; |
697 | gp_Pnt P; |
698 | gp_Vec DL1U, DL1V, DL2U, DL2V, DL2UV, DL3U, DL3UUV, DL3UVV, DL3V; |
699 | |
700 | if (theAlongU || theAlongV) |
701 | { |
702 | theMaxOrder = 0; |
703 | TColgp_Array2OfVec DerSurfL(0, theMaxOrder + theNU + 1, 0, theMaxOrder + theNV + 1); |
704 | switch (theMinOrder) |
705 | { |
706 | case 1: |
707 | theL->D1(theU, theV, P, DL1U, DL1V); |
708 | DerSurfL.SetValue(1, 0, DL1U); |
709 | DerSurfL.SetValue(0, 1, DL1V); |
710 | break; |
711 | case 2: |
712 | theL->D2(theU, theV, P, DL1U, DL1V, DL2U, DL2V, DL2UV); |
713 | DerSurfL.SetValue(1, 0, DL1U); |
714 | DerSurfL.SetValue(0, 1, DL1V); |
715 | DerSurfL.SetValue(1, 1, DL2UV); |
716 | DerSurfL.SetValue(2, 0, DL2U); |
717 | DerSurfL.SetValue(0, 2, DL2V); |
718 | break; |
719 | case 3: |
720 | theL->D3(theU, theV, P, DL1U, DL1V, DL2U, DL2V, DL2UV, DL3U, DL3V, DL3UUV, DL3UVV); |
721 | DerSurfL.SetValue(1, 0, DL1U); |
722 | DerSurfL.SetValue(0, 1, DL1V); |
723 | DerSurfL.SetValue(1, 1, DL2UV); |
724 | DerSurfL.SetValue(2, 0, DL2U); |
725 | DerSurfL.SetValue(0, 2, DL2V); |
726 | DerSurfL.SetValue(3, 0, DL3U); |
727 | DerSurfL.SetValue(2, 1, DL3UUV); |
728 | DerSurfL.SetValue(1, 2, DL3UVV); |
729 | DerSurfL.SetValue(0, 3, DL3V); |
730 | break; |
731 | default: |
732 | break; |
733 | } |
734 | |
735 | if (theNU <= theNV) |
736 | { |
737 | for (i = 0; i <= theMaxOrder + 1 + theNU; i++) |
738 | for (j = i; j <= theMaxOrder + theNV + 1; j++) |
739 | if (i + j > theMinOrder) |
740 | { |
741 | DerSurfL.SetValue(i, j, theL->DN(theU, theV, i, j)); |
742 | theDerSurf.SetValue(i, j, theBasisSurf->DN(theU, theV, i, j)); |
743 | if (i != j && j <= theNU + 1) |
744 | { |
745 | theDerSurf.SetValue(j, i, theBasisSurf->DN(theU, theV, j, i)); |
746 | DerSurfL.SetValue(j, i, theL->DN(theU, theV, j, i)); |
747 | } |
748 | } |
749 | } |
750 | else |
751 | { |
752 | for (j = 0; j <= theMaxOrder + 1 + theNV; j++) |
753 | for (i = j; i <= theMaxOrder + theNU + 1; i++) |
754 | if (i + j > theMinOrder) |
755 | { |
756 | DerSurfL.SetValue(i, j, theL->DN(theU, theV, i, j)); |
757 | theDerSurf.SetValue(i, j, theBasisSurf->DN(theU, theV, i, j)); |
758 | if (i != j && i <= theNV + 1) |
759 | { |
760 | theDerSurf.SetValue(j, i, theBasisSurf->DN(theU, theV, j, i)); |
761 | DerSurfL.SetValue(j, i, theL->DN(theU, theV, j, i)); |
762 | } |
763 | } |
764 | } |
765 | for (i = 0; i <= theMaxOrder + theNU; i++) |
766 | for (j = 0; j <= theMaxOrder + theNV; j++) |
767 | { |
768 | if (theAlongU) |
769 | theDerNUV.SetValue(i, j, CSLib::DNNUV(i, j, DerSurfL, theDerSurf)); |
770 | if (theAlongV) |
771 | theDerNUV.SetValue(i, j, CSLib::DNNUV(i, j, theDerSurf, DerSurfL)); |
772 | } |
773 | } |
774 | else |
775 | { |
776 | for (i = 0; i <= theMaxOrder + theNU+ 1; i++) |
777 | for (j = i; j <= theMaxOrder + theNV + 1; j++) |
778 | if (i + j > theMinOrder) |
779 | { |
780 | theDerSurf.SetValue(i, j, theBasisSurf->DN(theU, theV, i, j)); |
781 | if (i != j) |
782 | theDerSurf.SetValue(j, i, theBasisSurf->DN(theU, theV, j, i)); |
783 | } |
784 | for (i = 0; i <= theMaxOrder + theNU; i++) |
785 | for (j = 0; j <= theMaxOrder + theNV; j++) |
786 | theDerNUV.SetValue(i, j, CSLib::DNNUV(i, j, theDerSurf)); |
787 | } |
788 | } |
789 | |