src/GeomFill/GeomFill_LocationDraft.cxx

   1 // Created on: 1998-04-21
   2 // Created by: Stephanie HUMEAU
   3 // Copyright (c) 1998-1999 Matra Datavision
   4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
   5 //
   6 // This file is part of Open CASCADE Technology software library.
   7 //
   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.
  13 //
  14 // Alternatively, this file may be used under the terms of Open CASCADE
  15 // commercial license or contractual agreement.
  16
  17
  18 #include <Adaptor3d_HCurve.hxx>
  19 #include <Adaptor3d_HSurface.hxx>
  20 #include <Geom_Line.hxx>
  21 #include <Geom_Surface.hxx>
  22 #include <GeomAdaptor_HCurve.hxx>
  23 #include <GeomAdaptor_HSurface.hxx>
  24 #include <GeomFill_DraftTrihedron.hxx>
  25 #include <GeomFill_FunctionDraft.hxx>
  26 #include <GeomFill_LocationDraft.hxx>
  27 #include <GeomFill_LocationLaw.hxx>
  28 #include <GeomFill_Tensor.hxx>
  29 #include <GeomFill_TrihedronLaw.hxx>
  30 #include <gp_Dir.hxx>
  31 #include <gp_Mat.hxx>
  32 #include <gp_Pnt.hxx>
  33 #include <gp_Vec.hxx>
  34 #include <IntCurveSurface_HInter.hxx>
  35 #include <IntCurveSurface_Intersection.hxx>
  36 #include <IntCurveSurface_IntersectionPoint.hxx>
  37 #include <math_FunctionSetWithDerivatives.hxx>
  38 #include <math_Gauss.hxx>
  39 #include <math_Matrix.hxx>
  40 #include <math_NewtonFunctionSetRoot.hxx>
  41 #include <math_Vector.hxx>
  42 #include <Standard_NotImplemented.hxx>
  43 #include <Standard_OutOfRange.hxx>
  44 #include <Standard_Type.hxx>
  45
  46 IMPLEMENT_STANDARD_RTTIEXT(GeomFill_LocationDraft,GeomFill_LocationLaw)
  47
  48 //==================================================================
  49 //Function: GeomFill_LocationDraft
  50 //Purpose : constructor
  51 //==================================================================
  52 GeomFill_LocationDraft::GeomFill_LocationDraft
  53  (const gp_Dir& Direction,
  54   const Standard_Real Angle)
  55 {
  56   myDir = Direction; // direction de depouille
  57
  58   myAngle = Angle; // angle de depouille (teta prime)
  59
  60   mySurf.Nullify();
  61   myLaw = new (GeomFill_DraftTrihedron)(myDir, Angle); // triedre
  62   myNbPts = 41; // nb de points utilises pour les calculs
  63   myPoles2d = new (TColgp_HArray1OfPnt2d)(1, 2*myNbPts);
  64   Intersec = Standard_False; //intersection avec surface d'arret ?
  65   WithTrans = Standard_False;
  66 }
  67
  68
  69
  70 //==================================================================
  71 //Function: Copy
  72 //Purpose :
  73 //==================================================================
  74  Handle(GeomFill_LocationLaw) GeomFill_LocationDraft::Copy() const
  75 {
  76   Handle(GeomFill_TrihedronLaw) law;
  77   law = myLaw->Copy();
  78   Handle(GeomFill_LocationDraft) copy =
  79     new (GeomFill_LocationDraft) (myDir,myAngle);
  80   copy->SetCurve(myCurve);
  81   copy->SetStopSurf(mySurf);
  82   if (WithTrans) copy->SetTrsf(Trans);
  83
  84   return copy;
  85 }
  86
  87 //==================================================================
  88 //Function: SetTrsf
  89 //Purpose :
  90 //==================================================================
  91  void GeomFill_LocationDraft::SetTrsf(const gp_Mat& Transfo)
  92 {
  93   Trans = Transfo;
  94   gp_Mat Aux;
  95   Aux.SetIdentity();
  96   Aux -= Trans;
  97   WithTrans = Standard_False; // Au cas ou Trans = I
  98   for (Standard_Integer ii=1; ii<=3 && !WithTrans ; ii++)
  99     for (Standard_Integer jj=1; jj<=3 && !WithTrans; jj++)
 100       if (Abs(Aux.Value(ii, jj)) > 1.e-14)  WithTrans = Standard_True;
 101 }
 102
 103 //==================================================================
 104 //Function: SetCurve
 105 //Purpose : Calcul des poles sur la surfaces d'arret (intersection
 106 // entre la generatrice et la surface en myNbPts points de la section)
 107 //==================================================================
 108  void GeomFill_LocationDraft::SetCurve(const Handle(Adaptor3d_HCurve)& C)
 109 {
 110   myCurve = C;
 111   myTrimmed = C;
 112   myLaw->SetCurve(C);
 113
 114   Prepare();
 115 }
 116
 117 //==================================================================
 118 //Function: SetStopSurf
 119 //Purpose :
 120 //==================================================================
 121  void GeomFill_LocationDraft::SetStopSurf(const Handle(Adaptor3d_HSurface)& Surf)
 122 {
 123   mySurf = Surf;
 124   Prepare();
 125 }
 126
 127 //==================================================================
 128 //Function: SetAngle
 129 //Purpose :
 130 //==================================================================
 131  void GeomFill_LocationDraft::SetAngle(const Standard_Real Angle)
 132 {
 133   myAngle = Angle;
 134   myLaw->SetAngle(myAngle);
 135   Prepare();
 136 }
 137
 138 //==================================================================
 139 //Function: Prepare
 140 //Purpose : Poses les jalon de l'intersection : depouille / Surface
 141 //==================================================================
 142  void GeomFill_LocationDraft::Prepare()
 143 {
 144   if (mySurf.IsNull()) {
 145     Intersec = Standard_False;
 146     return;
 147   }
 148
 149   Intersec = Standard_True;
 150
 151   Standard_Integer ii,jj;
 152   Standard_Real f, l, t;
 153   gp_Pnt P;
 154   gp_Vec D,T,N,B;
 155   Handle(Geom_Line) L;
 156   IntCurveSurface_IntersectionPoint P1,P2;
 157   f = myCurve->FirstParameter();
 158   l = myCurve->LastParameter();
 159
 160   for (ii=1; ii<=myNbPts; ii++)
 161     {
 162       t = Standard_Real(myNbPts - ii)*f + Standard_Real(ii - 1)*l;
 163       t /= (myNbPts-1);
 164
 165       myCurve->D0(t, P);
 166       myLaw->D0(t,T,N,B);
 167
 168 // Generatrice
 169       D = Cos(myAngle)*B + Sin(myAngle)*N;
 170
 171       L = new (Geom_Line) (P, D);
 172
 173       IntCurveSurface_HInter Int; // intersection surface / generatrice
 174       Handle(GeomAdaptor_HCurve) AC = new (GeomAdaptor_HCurve) (L);
 175       Int.Perform(AC, mySurf); // calcul de l'intersection
 176
 177       if (Int.NbPoints() > 0) // il y a au moins 1 intersection
 178         {
 179           P1 = Int.Point(1);  // 1er  point d'intersection
 180
 181           for (jj=2 ; jj<=Int.NbPoints() ; jj++)
 182             {
 183               P2 = Int.Point(jj);
 184               if(P1.W() > P2.W()) P1 = P2; // point le plus proche
 185             }//for_jj
 186
 187           gp_Pnt2d p (P1.W(), t);        // point de la courbe
 188           gp_Pnt2d q (P1.U(),P1.V());    // point sur la surface
 189           myPoles2d->SetValue(2*ii-1,p); // point de la courbe (indice impair)
 190           myPoles2d->SetValue(2*ii,q);   // point sur la surface (indice pair)
 191         }
 192       else
 193         {// au moins un point ou il n'y a pas intersection
 194           Intersec = Standard_False;
 195         }
 196
 197     }//for_ii
 198 }
 199
 200 //==================================================================
 201 //Function: GetCurve
 202 //Purpose : return the path
 203 //==================================================================
 204  const Handle(Adaptor3d_HCurve)& GeomFill_LocationDraft::GetCurve() const
 205 {
 206   return myCurve;
 207 }
 208
 209
 210 //==================================================================
 211 //Function: D0
 212 //Purpose :
 213 //==================================================================
 214  Standard_Boolean GeomFill_LocationDraft::D0(const Standard_Real Param,
 215                                              gp_Mat& M,
 216                                              gp_Vec& V)
 217 {
 218   Standard_Boolean Ok;
 219   gp_Vec T,N,B;
 220   gp_Pnt P;
 221
 222   myTrimmed->D0(Param, P);
 223   V.SetXYZ(P.XYZ());
 224
 225   Ok = myLaw->D0(Param, T, N, B);
 226   if (!Ok) return Ok;
 227   M.SetCols(N.XYZ(), B.XYZ(), T.XYZ());
 228
 229   if (WithTrans) {
 230     M *= Trans;
 231   }
 232
 233   return Standard_True;
 234 }
 235
 236 //==================================================================
 237 //Function: D0
 238 //Purpose : calcul de l'intersection (C0) sur la surface
 239 //==================================================================
 240  Standard_Boolean GeomFill_LocationDraft::D0(const Standard_Real Param,
 241                                              gp_Mat& M,
 242                                              gp_Vec& V,
 243                                              TColgp_Array1OfPnt2d& Poles2d)
 244 {
 245   Standard_Boolean Ok;
 246 //  gp_Vec D,T,N,B,DT,DN,DB;
 247   gp_Vec D,T,N,B;
 248   gp_Pnt P;
 249
 250   myCurve->D0(Param, P);
 251   V.SetXYZ(P.XYZ());
 252   Ok = myLaw->D0(Param, T, N, B);
 253   if (!Ok) return Ok;
 254   M.SetCols(N.XYZ(), B.XYZ(), T.XYZ());
 255
 256   if (WithTrans) {
 257     M *= Trans;
 258   }
 259
 260   if (Intersec == Standard_True) {
 261     // la generatrice intersecte la surface d'arret
 262     // la generatrice
 263     D = Cos(myAngle)*B + Sin(myAngle)*N;
 264
 265     Handle(Geom_Line) L = new (Geom_Line) (P, D);
 266     Handle(GeomAdaptor_HCurve) G = new (GeomAdaptor_HCurve) (L);
 267
 268     Standard_Real t1,t2,Paramt1,t2Param;
 269     Standard_Real U0=0,V0=0,W0=0;
 270
 271     Standard_Integer ii = 1;
 272
 273     // on recherche l'intervalle auquel appartient Param
 274       while (ii<2*myNbPts && myPoles2d->Value(ii).Coord(2) < Param) ii=ii+2;
 275
 276     if (ii<2*myNbPts && !IsEqual(myPoles2d->Value(ii).Coord(2),Param))
 277       {
 278
 279         // interpolation lineaire pour initialiser le germe de la recherche
 280           t1 = myPoles2d->Value(ii).Coord(2);
 281           t2 = myPoles2d->Value(ii-2).Coord(2);
 282
 283           Paramt1 = (Param-t1) / (t2-t1);
 284           t2Param = (t2-Param) / (t2-t1);
 285
 286           W0 = myPoles2d->Value(ii-2).Coord(1)*Paramt1
 287                         + myPoles2d->Value(ii).Coord(1)*t2Param;
 288           U0 = myPoles2d->Value(ii-1).Coord(1)*Paramt1
 289                         + myPoles2d->Value(ii+1).Coord(1)*t2Param;
 290           V0 = myPoles2d->Value(ii-1).Coord(2)*Paramt1
 291                         + myPoles2d->Value(ii+1).Coord(2)*t2Param;
 292         }//if
 293       // on est sur un param ou les points ont deja ete calcules
 294       else if (ii<2*myNbPts && IsEqual(myPoles2d->Value(ii).Coord(2) ,Param))
 295         {
 296           W0 = myPoles2d->Value(ii).Coord(1);
 297           U0 = myPoles2d->Value(ii+1).Coord(1);
 298           V0 = myPoles2d->Value(ii+1).Coord(2);
 299         }//else if
 300
 301  // recherche de la solution (pt d'intersection generatrice / surface)
 302       // point initial
 303       math_Vector X(1,3);
 304       X(1) = W0;
 305       X(2) = U0;
 306       X(3) = V0;
 307
 308       // tolerance sur X
 309       math_Vector XTol(1,3);
 310       XTol.Init(0.00001);
 311
 312       // tolerance sur F
 313       Standard_Real FTol = 0.0000001;
 314       Standard_Integer Iter = 100;
 315
 316       // fonction dont il faut trouver la racine : G(W)-S(U,V)=0
 317       GeomFill_FunctionDraft E(mySurf , G);
 318
 319       // resolution
 320       math_NewtonFunctionSetRoot Result(E, XTol, FTol, Iter);
 321       Result.Perform(E, X);
 322
 323       if (Result.IsDone())
 324       {
 325         math_Vector R(1,3);
 326         Result.Root(R);    // solution
 327
 328         gp_Pnt2d p (R(2), R(3));  // point sur la surface
 329         gp_Pnt2d q (R(1), Param); // point de la courbe
 330         Poles2d.SetValue(1,p);
 331         Poles2d.SetValue(2,q);
 332       }
 333       else {
 334         return Standard_False;
 335       }
 336   }// if_Intersec
 337
 338   // la generatrice n'intersecte pas la surface d'arret
 339   return Standard_True;
 340  }
 341
 342 //==================================================================
 343 //Function: D1
 344 //Purpose : calcul de l'intersection (C1) sur la surface
 345 //==================================================================
 346  Standard_Boolean GeomFill_LocationDraft::D1(const Standard_Real Param,
 347                                              gp_Mat& M,
 348                                              gp_Vec& V,
 349                                              gp_Mat& DM,
 350                                              gp_Vec& DV,
 351                                              TColgp_Array1OfPnt2d& Poles2d,
 352                                              TColgp_Array1OfVec2d& DPoles2d)
 353 {
 354   Standard_Boolean Ok;
 355   gp_Vec D,T,N,B,DT,DN,DB;
 356   gp_Pnt P;
 357
 358   myCurve->D1(Param, P, DV);
 359   V.SetXYZ(P.XYZ());
 360
 361   Ok = myLaw->D1(Param, T, DT, N, DN, B, DB);
 362   if (!Ok) return Standard_False;
 363
 364   M.SetCols(N.XYZ(), B.XYZ(), T.XYZ());
 365   DM.SetCols(DN.XYZ(), DB.XYZ(), DT.XYZ());
 366
 367   if (WithTrans) {
 368     M *= Trans;
 369     DM *= Trans;
 370   }
 371
 372   if (Intersec == Standard_True)
 373     { // la generatrice intersecte la surface d'arret
 374       // la generatrice
 375       D = Cos(myAngle)*B + Sin(myAngle)*N;
 376
 377       Handle(Geom_Line) L = new (Geom_Line) (P, D);
 378       Handle(GeomAdaptor_HCurve) G = new (GeomAdaptor_HCurve) (L);
 379
 380       Standard_Real t1,t2,Paramt1,t2Param;
 381       Standard_Real U0=0,V0=0,W0=0;
 382
 383       Standard_Integer ii = 1;
 384
 385       // recherche de la solution (pt d'intersection  generatrice / surface)
 386
 387       // on recherche l'intervalle auquel appartient Param
 388       while (ii < 2*myNbPts && myPoles2d->Value(ii).Coord(2) < Param) ii=ii+2;
 389
 390
 391       if (ii < 2*myNbPts && !IsEqual(myPoles2d->Value(ii).Coord(2) ,Param))
 392         {
 393           // interpolation lineaire pour initialiser le germe de la recherche
 394           t1 = myPoles2d->Value(ii).Coord(2);
 395           t2 = myPoles2d->Value(ii-2).Coord(2);
 396
 397           Paramt1 = (Param-t1) / (t2-t1);
 398           t2Param = (t2-Param) / (t2-t1);
 399
 400           W0 = myPoles2d->Value(ii-2).Coord(1)*Paramt1
 401                        + myPoles2d->Value(ii).Coord(1)*t2Param;
 402           U0 = myPoles2d->Value(ii-1).Coord(1)*Paramt1
 403                        + myPoles2d->Value(ii+1).Coord(1)*t2Param;
 404           V0 = myPoles2d->Value(ii-1).Coord(2)*Paramt1
 405                        + myPoles2d->Value(ii+1).Coord(2)*t2Param;
 406         }//if
 407       else if (ii<2*myNbPts && IsEqual(myPoles2d->Value(ii).Coord(2) ,Param))
 408         {
 409           W0 = myPoles2d->Value(ii).Coord(1);
 410           U0 = myPoles2d->Value(ii+1).Coord(1);
 411           V0 = myPoles2d->Value(ii+1).Coord(2);
 412         }//else if
 413
 414       // germe
 415       math_Vector X(1,3);
 416       X(1) = W0;
 417       X(2) = U0;
 418       X(3) = V0;
 419
 420       // tolerance sur X
 421       math_Vector XTol(1,3);
 422       XTol.Init(0.0001);
 423
 424       // tolerance sur F
 425       Standard_Real FTol = 0.000001;
 426       Standard_Integer Iter = 100;
 427
 428
 429       // fonction dont il faut trouver la racine : G(W)-S(U,V)=0
 430       GeomFill_FunctionDraft E(mySurf,G);
 431
 432
 433       // resolution
 434       math_NewtonFunctionSetRoot Result(E, XTol, FTol, Iter);
 435       Result.Perform(E, X);
 436
 437       if (Result.IsDone())
 438       {
 439         math_Vector R(1,3);
 440         Result.Root(R);    // solution
 441
 442         gp_Pnt2d p (R(2), R(3));  // point sur la surface
 443         gp_Pnt2d q (R(1), Param); // point de la courbe
 444         Poles2d.SetValue(1,p);
 445         Poles2d.SetValue(2,q);
 446
 447         // derivee de la fonction par rapport a Param
 448         math_Vector DEDT(1,3,0);
 449         E.DerivT(myTrimmed, Param, R(1), DN, myAngle, DEDT); // dE/dt => DEDT
 450
 451         math_Vector DSDT (1,3,0);
 452         math_Matrix DEDX (1,3,1,3,0);
 453         E.Derivatives(R, DEDX);  // dE/dx au point R => DEDX
 454
 455         // resolution du syst. lin. : DEDX*DSDT = -DEDT
 456         math_Gauss Ga(DEDX);
 457         if (Ga.IsDone())
 458         {
 459           Ga.Solve (DEDT.Opposite(), DSDT); // resolution du syst. lin.
 460           gp_Vec2d dp (DSDT(2), DSDT(3));    // surface
 461           gp_Vec2d dq (DSDT(1), 1);          //  courbe
 462           DPoles2d.SetValue(1, dp);
 463           DPoles2d.SetValue(2, dq);
 464         }//if
 465
 466
 467       }//if_Result
 468       else {// la generatrice n'intersecte pas la surface d'arret
 469         return Standard_False;
 470       }
 471   }// if_Intersec
 472   return Standard_True;
 473 }
 474
 475 //==================================================================
 476 //Function: D2
 477 //Purpose : calcul de l'intersection (C2) sur la surface
 478 //==================================================================
 479  Standard_Boolean GeomFill_LocationDraft::D2(const Standard_Real Param,
 480                                              gp_Mat& M,
 481                                              gp_Vec& V,
 482                                              gp_Mat& DM,
 483                                              gp_Vec& DV,
 484                                              gp_Mat& D2M,
 485                                              gp_Vec& D2V,
 486                                              TColgp_Array1OfPnt2d& Poles2d,
 487                                              TColgp_Array1OfVec2d& DPoles2d,
 488                                              TColgp_Array1OfVec2d& D2Poles2d)
 489 {
 490   Standard_Boolean Ok;
 491   gp_Vec D,T,N,B,DT,DN,DB,D2T,D2N,D2B;
 492   gp_Pnt P;
 493
 494   myCurve->D2(Param, P, DV, D2V);
 495   V.SetXYZ(P.XYZ());
 496
 497   Ok = myLaw->D2(Param, T, DT, D2T, N, DN, D2N, B, DB, D2B);
 498   if (!Ok) return Ok;
 499
 500   M.SetCols(N.XYZ(), B.XYZ(), T.XYZ());
 501   DM.SetCols(DN.XYZ(), DB.XYZ(), DT.XYZ());
 502   D2M.SetCols(D2N.XYZ(), D2B.XYZ(), D2T.XYZ());
 503
 504   if (WithTrans) {
 505     M *= Trans;
 506     DM *= Trans;
 507     D2M *= Trans;
 508   }
 509   if (Intersec == Standard_True)
 510     {// la generatrice intersecte la surface d'arret
 511
 512       // la generatrice
 513       D = Cos(myAngle) * B + Sin(myAngle) * N;
 514
 515       Handle(Geom_Line) L = new (Geom_Line) (P, D);
 516       Handle(GeomAdaptor_HCurve) G = new (GeomAdaptor_HCurve) (L);
 517
 518       Standard_Real t1,t2,Paramt1,t2Param;
 519       Standard_Real U0=0,V0=0,W0=0;
 520
 521       Standard_Integer ii = 1;
 522
 523       // on recherche l'intervalle auquel appartient Param
 524       while (ii<2*myNbPts && myPoles2d->Value(ii).Coord(2) < Param) ii=ii+2;
 525
 526       if (ii<2*myNbPts && !IsEqual(myPoles2d->Value(ii).Coord(2) ,Param))
 527         {
 528
 529           // interpolation lineaire pour initialiser le germe de la recherche
 530           t1 = myPoles2d->Value(ii).Coord(2);
 531           t2 = myPoles2d->Value(ii-2).Coord(2);
 532
 533           Paramt1 = (Param-t1) / (t2-t1);
 534           t2Param = (t2-Param) / (t2-t1);
 535
 536           W0 = myPoles2d->Value(ii-2).Coord(1)*Paramt1 +
 537             myPoles2d->Value(ii).Coord(1)*t2Param;
 538           U0 = myPoles2d->Value(ii-1).Coord(1)*Paramt1 +
 539             myPoles2d->Value(ii+1).Coord(1)*t2Param;
 540           V0 = myPoles2d->Value(ii-1).Coord(2)*Paramt1 +
 541             myPoles2d->Value(ii+1).Coord(2)*t2Param;
 542         }//if
 543       else if (ii<2*myNbPts  && IsEqual(myPoles2d->Value(ii).Coord(2) ,Param))
 544         {
 545           W0 = myPoles2d->Value(ii).Coord(1);
 546           U0 = myPoles2d->Value(ii+1).Coord(1);
 547           V0 = myPoles2d->Value(ii+1).Coord(2);
 548         }//else if
 549
 550 // recherche de la solution (pt d'intersection generatrice / surface)
 551       // germe
 552       math_Vector X(1,3);
 553       X(1) = W0;
 554       X(2) = U0;
 555       X(3) = V0;
 556
 557       // tolerance sur X
 558       math_Vector XTol(1,3);
 559       XTol.Init(0.0001);
 560
 561       // tolerance sur F
 562       Standard_Real FTol = 0.000001;
 563       Standard_Integer Iter = 150;
 564
 565
 566       // fonction dont il faut trouver la racine : G(W)-S(U,V)=0
 567       GeomFill_FunctionDraft E(mySurf,G);
 568
 569
 570       // resolution
 571       math_NewtonFunctionSetRoot Result (E, XTol, FTol, Iter);
 572       Result.Perform(E, X);
 573
 574       if (Result.IsDone())
 575         {
 576           math_Vector R(1,3);
 577           Result.Root(R);    // solution
 578
 579           // solution
 580           gp_Pnt2d p (R(2), R(3));  // point sur la surface
 581           gp_Pnt2d q (R(1), Param); // point de la courbe
 582           Poles2d.SetValue(1,p);
 583           Poles2d.SetValue(2,q);
 584
 585
 586      // premiere derivee de la fonction
 587           math_Vector DEDT(1,3,0);
 588           E.DerivT(myTrimmed, Param, R(1), DN, myAngle, DEDT); // dE/dt => DEDT
 589
 590           math_Vector DSDT (1,3,0);
 591           math_Matrix DEDX (1,3,1,3,0);
 592           E.Derivatives(R, DEDX);  // dE/dx => DEDX
 593
 594           // resolution du syst. lin.
 595           math_Gauss Ga (DEDX);
 596           if (Ga.IsDone())
 597             {
 598               Ga.Solve (DEDT.Opposite(), DSDT);
 599               gp_Vec2d dp (DSDT(2), DSDT(3));   // surface
 600               gp_Vec2d dq (DSDT(1), 1);         //  courbe
 601               DPoles2d.SetValue(1, dp);
 602               DPoles2d.SetValue(2, dq);
 603             }//if
 604
 605
 606      // deuxieme derivee
 607           GeomFill_Tensor D2EDX2(3,3,3);
 608           E.Deriv2X(R, D2EDX2); // d2E/dx2
 609
 610           math_Vector D2EDT2(1,3,0);
 611           E.Deriv2T(myTrimmed, Param, R(1), D2N, myAngle ,D2EDT2); // d2E/dt2
 612
 613           math_Matrix D2EDTDX(1,3,1,3,0);
 614           E.DerivTX(DN, myAngle, D2EDTDX); // d2E/dtdx
 615
 616           math_Vector D2SDT2(1,3,0); // d2s/dt2
 617           math_Matrix aT(1,3,1,3,0);
 618           D2EDX2.Multiply(DSDT,aT);
 619
 620           // resolution du syst. lin.
 621           math_Gauss Ga1 (DEDX);
 622           if (Ga1.IsDone())
 623             {
 624               Ga1.Solve ( -aT*DSDT - 2*D2EDTDX*DSDT - D2EDT2 , D2SDT2);
 625               gp_Vec2d d2p (D2SDT2(2), D2SDT2(3));  // surface
 626               gp_Vec2d d2q (D2SDT2(1), 0);          // courbe
 627               D2Poles2d.SetValue(1, d2p);
 628               D2Poles2d.SetValue(2, d2q);
 629             }//if
 630           else {// la generatrice n'intersecte pas la surface d'arret
 631             return Standard_False;
 632           }
 633         }//if_Result
 634     } //if_Intersec
 635
 636   return Standard_True;
 637 }
 638
 639 //==================================================================
 640 //Function : HasFirstRestriction
 641 //Purpose :
 642 //==================================================================
 643  Standard_Boolean GeomFill_LocationDraft::HasFirstRestriction() const
 644 {
 645   return Standard_False;
 646 }
 647
 648 //==================================================================
 649 //Function : HasLastRestriction
 650 //Purpose :
 651 //==================================================================
 652  Standard_Boolean GeomFill_LocationDraft::HasLastRestriction() const
 653 {
 654
 655   if (Intersec == Standard_True) return Standard_True;
 656   else return Standard_False;
 657 }
 658
 659 //==================================================================
 660 //Function : TraceNumber
 661 //Purpose :
 662 //==================================================================
 663  Standard_Integer GeomFill_LocationDraft::TraceNumber() const
 664 {
 665   if (Intersec == Standard_True) return 1;
 666   else return 0;
 667 }
 668
 669 //==================================================================
 670 //Function:NbIntervals
 671 //Purpose :
 672 //==================================================================
 673  Standard_Integer GeomFill_LocationDraft::NbIntervals
 674  (const GeomAbs_Shape S) const
 675 {
 676   return myLaw->NbIntervals(S);
 677 }
 678
 679 //==================================================================
 680 //Function:Intervals
 681 //Purpose :
 682 //==================================================================
 683  void GeomFill_LocationDraft::Intervals(TColStd_Array1OfReal& T,
 684                                         const GeomAbs_Shape S) const
 685 {
 686   myLaw->Intervals(T, S);
 687 }
 688
 689 //==================================================================
 690 //Function:SetInterval
 691 //Purpose :
 692 //==================================================================
 693  void GeomFill_LocationDraft::SetInterval(const Standard_Real First,
 694                                           const Standard_Real Last)
 695 {
 696   myLaw->SetInterval( First, Last);
 697   myTrimmed = myCurve->Trim( First, Last, 0);
 698 }
 699 //==================================================================
 700 //Function: GetInterval
 701 //Purpose :
 702 //==================================================================
 703  void GeomFill_LocationDraft::GetInterval(Standard_Real& First,
 704                                           Standard_Real& Last) const
 705 {
 706   First = myTrimmed->FirstParameter();
 707   Last = myTrimmed->LastParameter();
 708 }
 709
 710 //==================================================================
 711 //Function: GetDomain
 712 //Purpose :
 713 //==================================================================
 714  void GeomFill_LocationDraft::GetDomain(Standard_Real& First,
 715                                         Standard_Real& Last) const
 716 {
 717   First = myCurve->FirstParameter();
 718   Last = myCurve->LastParameter();
 719 }
 720
 721 //==================================================================
 722 //function : Resolution
 723 //purpose  :
 724 //==================================================================
 725 void GeomFill_LocationDraft::Resolution (const Standard_Integer Index,
 726                                          const Standard_Real Tol,
 727                                          Standard_Real& TolU,
 728                                          Standard_Real& TolV) const
 729
 730 {
 731   if (Index==1) {
 732       TolU = mySurf->UResolution(Tol);
 733       TolV = mySurf->VResolution(Tol);
 734     }
 735   else {
 736       TolU = Tol;
 737       TolV = Tol;
 738     }
 739 }
 740
 741 //==================================================================
 742 //Function:GetMaximalNorm
 743 //Purpose :  On suppose les triedres normes => return 1
 744 //==================================================================
 745  Standard_Real GeomFill_LocationDraft::GetMaximalNorm()
 746 {
 747   return 1.;
 748 }
 749
 750 //==================================================================
 751 //Function:GetAverageLaw
 752 //Purpose :
 753 //==================================================================
 754  void GeomFill_LocationDraft::GetAverageLaw(gp_Mat& AM,
 755                                             gp_Vec& AV)
 756 {
 757   Standard_Integer ii;
 758   Standard_Real U, delta;
 759   gp_Vec V1,V2,V3, V;
 760
 761   myLaw->GetAverageLaw(V1, V2, V3);
 762   AM.SetCols(V1.XYZ(), V2.XYZ(), V3.XYZ());
 763
 764   AV.SetCoord(0., 0., 0.);
 765   delta = (myTrimmed->LastParameter() - myTrimmed->FirstParameter())/10;
 766   U=  myTrimmed->FirstParameter();
 767   for (ii=0; ii<=10; ii++, U+=delta) {
 768     V.SetXYZ( myTrimmed->Value(U).XYZ() );
 769     AV += V;
 770   }
 771   AV /= 11;
 772 }
 773
 774 //==================================================================
 775 //Function : IsTranslation
 776 //Purpose :
 777 //==================================================================
 778 // Standard_Boolean GeomFill_LocationDraft::IsTranslation(Standard_Real& Error) const
 779  Standard_Boolean GeomFill_LocationDraft::IsTranslation(Standard_Real& ) const
 780 {
 781   return myLaw->IsConstant();
 782 }
 783
 784 //==================================================================
 785 //Function : IsRotation
 786 //Purpose :
 787 //==================================================================
 788  Standard_Boolean GeomFill_LocationDraft::IsRotation(Standard_Real& Error)  const
 789 {
 790   GeomAbs_CurveType Type;
 791   Error = 0;
 792   Type = myCurve->GetType();
 793   if (Type == GeomAbs_Circle) {
 794     return myLaw->IsOnlyBy3dCurve();
 795   }
 796   return Standard_False;
 797 }
 798
 799 //==================================================================
 800 //Function : Rotation
 801 //Purpose :
 802 //==================================================================
 803  void GeomFill_LocationDraft::Rotation(gp_Pnt& Centre)  const
 804 {
 805   Centre =  myCurve->Circle().Location();
 806 }
 807
 808
 809 //==================================================================
 810 //Function : IsIntersec
 811 //Purpose :
 812 //==================================================================
 813  Standard_Boolean GeomFill_LocationDraft::IsIntersec()  const
 814 {
 815  return Intersec;
 816 }
 817
 818
 819 //==================================================================
 820 //Function : Direction
 821 //Purpose :
 822 //==================================================================
 823  gp_Dir GeomFill_LocationDraft::Direction()  const
 824 {
 825  return myDir;
 826 }