1 // Created on: 1993-12-21
2 // Created by: Jacques GOUSSARD
3 // Copyright (c) 1993-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
6 // This file is part of Open CASCADE Technology software library.
8 // This library is free software; you can redistribute it and/or modify it under
9 // the terms of the GNU Lesser General Public License version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
18 #include <Adaptor2d_HCurve2d.hxx>
19 #include <Adaptor3d_HCurve.hxx>
20 #include <Adaptor3d_HSurface.hxx>
21 #include <BlendFunc.hxx>
22 #include <BlendFunc_EvolRadInv.hxx>
23 #include <Law_Function.hxx>
24 #include <math_Matrix.hxx>
25 #include <Precision.hxx>
30 BlendFunc_EvolRadInv::BlendFunc_EvolRadInv(const Handle(Adaptor3d_HSurface)& S1,
31 const Handle(Adaptor3d_HSurface)& S2,
32 const Handle(Adaptor3d_HCurve)& C,
33 const Handle(Law_Function)& Law) :
34 surf1(S1),surf2(S2),curv(C)
39 void BlendFunc_EvolRadInv::Set(const Standard_Integer Choix)
77 void BlendFunc_EvolRadInv::Set(const Standard_Boolean OnFirst,
78 const Handle(Adaptor2d_HCurve2d)& C)
84 Standard_Integer BlendFunc_EvolRadInv::NbEquations () const
90 void BlendFunc_EvolRadInv::GetTolerance(math_Vector& Tolerance,
91 const Standard_Real Tol) const
93 Tolerance(1) = csurf->Resolution(Tol);
94 Tolerance(2) = curv->Resolution(Tol);
96 Tolerance(3) = surf2->UResolution(Tol);
97 Tolerance(4) = surf2->VResolution(Tol);
100 Tolerance(3) = surf1->UResolution(Tol);
101 Tolerance(4) = surf1->VResolution(Tol);
106 void BlendFunc_EvolRadInv::GetBounds(math_Vector& InfBound,
107 math_Vector& SupBound) const
109 InfBound(1) = csurf->FirstParameter();
110 InfBound(2) = curv->FirstParameter();
111 SupBound(1) = csurf->LastParameter();
112 SupBound(2) = curv->LastParameter();
115 InfBound(3) = surf2->FirstUParameter();
116 InfBound(4) = surf2->FirstVParameter();
117 SupBound(3) = surf2->LastUParameter();
118 SupBound(4) = surf2->LastVParameter();
119 if(!Precision::IsInfinite(InfBound(3)) &&
120 !Precision::IsInfinite(SupBound(3))) {
121 const Standard_Real range = (SupBound(3) - InfBound(3));
122 InfBound(3) -= range;
123 SupBound(3) += range;
125 if(!Precision::IsInfinite(InfBound(4)) &&
126 !Precision::IsInfinite(SupBound(4))) {
127 const Standard_Real range = (SupBound(4) - InfBound(4));
128 InfBound(4) -= range;
129 SupBound(4) += range;
133 InfBound(3) = surf1->FirstUParameter();
134 InfBound(4) = surf1->FirstVParameter();
135 SupBound(3) = surf1->LastUParameter();
136 SupBound(4) = surf1->LastVParameter();
137 if(!Precision::IsInfinite(InfBound(3)) &&
138 !Precision::IsInfinite(SupBound(3))) {
139 const Standard_Real range = (SupBound(3) - InfBound(3));
140 InfBound(3) -= range;
141 SupBound(3) += range;
143 if(!Precision::IsInfinite(InfBound(4)) &&
144 !Precision::IsInfinite(SupBound(4))) {
145 const Standard_Real range = (SupBound(4) - InfBound(4));
146 InfBound(4) -= range;
147 SupBound(4) += range;
152 Standard_Boolean BlendFunc_EvolRadInv::IsSolution(const math_Vector& Sol,
153 const Standard_Real Tol)
155 math_Vector valsol(1,4);
157 if (Abs(valsol(1)) <= Tol &&
158 (valsol(2)*valsol(2) + valsol(3)*valsol(3) + valsol(4)*valsol(4)) <= Tol*Tol)
159 return Standard_True;
161 return Standard_False;
165 Standard_Boolean BlendFunc_EvolRadInv::Value(const math_Vector& X,
168 const Standard_Real ray = fevol->Value(X(2));
172 curv->D1(X(2),ptcur,d1cur);
174 const gp_Vec nplan = d1cur.Normalized();
175 const Standard_Real theD = -(nplan.XYZ().Dot(ptcur.XYZ()));
177 const gp_Pnt2d pt2d(csurf->Value(X(1)));
180 gp_Vec d1u1,d1v1,d1u2,d1v2;
183 surf1->D1(pt2d.X(),pt2d.Y(),pts1,d1u1,d1v1);
184 surf2->D1(X(3),X(4),pts2,d1u2,d1v2);
188 surf1->D1(X(3),X(4),pts1,d1u1,d1v1);
189 surf2->D1(pt2d.X(),pt2d.Y(),pts2,d1u2,d1v2);
192 F(1) = (nplan.X() * (pts1.X() + pts2.X()) +
193 nplan.Y() * (pts1.Y() + pts2.Y()) +
194 nplan.Z() * (pts1.Z() + pts2.Z())) /2. + theD;
196 gp_Vec ns1 = d1u1.Crossed(d1v1);
197 if (ns1.Magnitude() < Eps) {
198 if (first) BlendFunc::ComputeNormal(surf1, pt2d, ns1);
200 gp_Pnt2d P(X(3), X(4));
201 BlendFunc::ComputeNormal(surf1, P, ns1);
205 gp_Vec ns2 = d1u2.Crossed(d1v2).XYZ();
206 if (ns2.Magnitude() < Eps) {
207 if (!first) BlendFunc::ComputeNormal(surf2, pt2d, ns2);
209 gp_Pnt2d P(X(3), X(4));
210 BlendFunc::ComputeNormal(surf2, P, ns2);
214 const gp_Vec ncrossns1 = nplan.Crossed(ns1);
215 const gp_Vec ncrossns2 = nplan.Crossed(ns2);
216 Standard_Real norm1 = ncrossns1.Magnitude();
217 Standard_Real norm2 = ncrossns2.Magnitude();
227 const Standard_Real ndotns1 = nplan.Dot(ns1);
228 const Standard_Real ndotns2 = nplan.Dot(ns2);
229 ns1.SetLinearForm(ndotns1/norm1,nplan, -1./norm1,ns1);
230 ns2.SetLinearForm(ndotns2/norm2,nplan, -1./norm2,ns2);
231 resul.SetLinearForm(sg1*ray,ns1,-1.,pts2.XYZ(),-sg2*ray,ns2,pts1.XYZ());
236 return Standard_True;
239 Standard_Boolean BlendFunc_EvolRadInv::Derivatives(const math_Vector& X,
243 return Values (X, F, D);
246 Standard_Boolean BlendFunc_EvolRadInv::Values(const math_Vector& X,
250 Standard_Real ray,dray;
251 fevol->D1(X(2),ray,dray);
255 curv->D2(X(2),ptcur,d1cur,d2cur);
257 const Standard_Real normtgcur = d1cur.Magnitude();
258 const gp_Vec nplan = d1cur.Normalized();
259 const Standard_Real theD = -(nplan.XYZ().Dot(ptcur.XYZ()));
262 dnplan.SetLinearForm(-nplan.Dot(d2cur),nplan,d2cur);
267 csurf->D1(X(1),p2d,v2d);
270 gp_Vec d1u1,d1v1,d1u2,d1v2;
271 gp_Vec d2u1,d2v1,d2u2,d2v2,d2uv1,d2uv2;
275 surf1->D2(p2d.X(),p2d.Y(),pts1,d1u1,d1v1,d2u1,d2v1,d2uv1);
276 surf2->D2(X(3),X(4),pts2,d1u2,d1v2,d2u2,d2v2,d2uv2);
278 temp.SetLinearForm(v2d.X(),d1u1,v2d.Y(),d1v1);
279 D(1,1) = nplan.Dot(temp)/2.;
280 temp.SetXYZ(0.5*(pts1.XYZ()+pts2.XYZ())-ptcur.XYZ());
281 D(1,2) = dnplan.Dot(temp) - normtgcur;
282 D(1,3) = nplan.Dot(d1u2)/2.;
283 D(1,4) = nplan.Dot(d1v2)/2.;
287 surf1->D2(X(3),X(4),pts1,d1u1,d1v1,d2u1,d2v1,d2uv1);
288 surf2->D2(p2d.X(),p2d.Y(),pts2,d1u2,d1v2,d2u2,d2v2,d2uv2);
290 temp.SetLinearForm(v2d.X(),d1u2,v2d.Y(),d1v2);
291 D(1,1) = nplan.Dot(temp)/2.;
292 temp.SetXYZ(0.5*(pts1.XYZ()+pts2.XYZ())-ptcur.XYZ());
293 D(1,2) = dnplan.Dot(temp) - normtgcur;
294 D(1,3) = nplan.Dot(d1u1)/2.;
295 D(1,4) = nplan.Dot(d1v1)/2.;
298 F(1) = (nplan.X()* (pts1.X() + pts2.X()) +
299 nplan.Y()* (pts1.Y() + pts2.Y()) +
300 nplan.Z()* (pts1.Z() + pts2.Z())) /2. + theD;
302 gp_Vec ns1 = d1u1.Crossed(d1v1);
303 if (ns1.Magnitude() < Eps) {
304 if (first) BlendFunc::ComputeNormal(surf1, p2d, ns1);
306 gp_Pnt2d P(X(3), X(4));
307 BlendFunc::ComputeNormal(surf1, P, ns1);
311 gp_Vec ns2 = d1u2.Crossed(d1v2);
312 if (ns2.Magnitude() < Eps) {
313 if (!first) BlendFunc::ComputeNormal(surf2, p2d, ns2);
315 gp_Pnt2d P(X(3), X(4));
316 BlendFunc::ComputeNormal(surf2, P, ns2);
320 const gp_Vec ncrossns1 = nplan.Crossed(ns1);
321 const gp_Vec ncrossns2 = nplan.Crossed(ns2);
322 Standard_Real norm1 = ncrossns1.Magnitude();
323 Standard_Real norm2 = ncrossns2.Magnitude();
331 gp_Vec resul1,resul2;
332 Standard_Real grosterme;
334 const Standard_Real ndotns1 = nplan.Dot(ns1);
335 const Standard_Real ndotns2 = nplan.Dot(ns2);
336 temp.SetLinearForm(ndotns1/norm1,nplan, -1./norm1,ns1);
337 resul1.SetLinearForm(sg1*ray,temp,gp_Vec(pts2,pts1));
338 temp.SetLinearForm(ndotns2/norm2,nplan,-1./norm2,ns2);
339 resul1.Subtract(sg2*ray*temp);
345 // Derivee par rapport a u1
347 temp = d2u1.Crossed(d1v1).Added(d1u1.Crossed(d2uv1));
348 grosterme = ncrossns1.Dot(nplan.Crossed(temp))/norm1/norm1;
350 (-sg1*ray/norm1*(grosterme*ndotns1-nplan.Dot(temp)),nplan,
351 sg1*ray*grosterme/norm1,ns1,
356 // Derivee par rapport a v1
358 temp = d2uv1.Crossed(d1v1).Added(d1u1.Crossed(d2v1));
359 grosterme = ncrossns1.Dot(nplan.Crossed(temp))/norm1/norm1;
361 (-sg1*ray/norm1*(grosterme*ndotns1-nplan.Dot(temp)),nplan,
362 sg1*ray*grosterme/norm1,ns1,
367 D(2,1) = v2d.X()*resul1.X() + v2d.Y()*resul2.X();
368 D(3,1) = v2d.X()*resul1.Y() + v2d.Y()*resul2.Y();
369 D(4,1) = v2d.X()*resul1.Z() + v2d.Y()*resul2.Z();
381 // derivee par rapport a w (parametre sur ligne guide)
383 grosterme = ncrossns1.Dot(dnplan.Crossed(ns1))/norm1/norm1;
384 resul1.SetLinearForm(-sg1/norm1*(grosterme*ndotns1-dnplan.Dot(ns1)),nplan,
385 sg1*ndotns1/norm1,dnplan,
386 sg1*grosterme/norm1,ns1);
389 grosterme = ncrossns2.Dot(dnplan.Crossed(ns2))/norm2/norm2;
390 resul2.SetLinearForm(sg2/norm2*(grosterme*ndotns2-dnplan.Dot(ns2)),nplan,
391 -sg2*ndotns2/norm2,dnplan,
392 -sg2*grosterme/norm2,ns2);
395 D(2,2) = ray*(resul1.X() + resul2.X());
396 D(3,2) = ray*(resul1.Y() + resul2.Y());
397 D(4,2) = ray*(resul1.Z() + resul2.Z());
399 temp.SetLinearForm(sg1*ndotns1/norm1 - sg2*ndotns2/norm2,nplan,
403 D(2,2) += dray*temp.X();
404 D(3,2) += dray*temp.Y();
405 D(4,2) += dray*temp.Z();
409 // Derivee par rapport a u2
410 temp = d2u2.Crossed(d1v2).Added(d1u2.Crossed(d2uv2));
411 grosterme = ncrossns2.Dot(nplan.Crossed(temp))/norm2/norm2;
412 resul1.SetLinearForm(sg2*ray/norm2*(grosterme*ndotns2-nplan.Dot(temp)),nplan,
413 -sg2*ray*grosterme/norm2,ns2,
415 resul1.Subtract(d1u2);
417 // Derivee par rapport a v2
418 temp = d2uv2.Crossed(d1v2).Added(d1u2.Crossed(d2v2));
419 grosterme = ncrossns2.Dot(nplan.Crossed(temp))/norm2/norm2;
420 resul2.SetLinearForm(sg2*ray/norm2*(grosterme*ndotns2-nplan.Dot(temp)),nplan,
421 -sg2*ray*grosterme/norm2,ns2,
423 resul2.Subtract(d1v2);
426 D(2,1) = v2d.X()*resul1.X() + v2d.Y()*resul2.X();
427 D(3,1) = v2d.X()*resul1.Y() + v2d.Y()*resul2.Y();
428 D(4,1) = v2d.X()*resul1.Z() + v2d.Y()*resul2.Z();
440 return Standard_True;