0022781: Regression in GCPnts_TangentialDeflection
[occt.git] / src / GeometryTest / GeometryTest_CurveCommands.cxx
1 // File:        DrawTrSurf_1.cxx
2 // Created:     Thu Aug 12 19:33:03 1993
3 // Author:      Bruno DUMORTIER
4 //              <dub@topsn3>
5 // 09/06/97 : JPI : suppression des commandes redondantes suite a la creation de GeomliteTest
6
7 #include <GeometryTest.hxx>
8 #include <Draw_Appli.hxx>
9 #include <DrawTrSurf.hxx>
10 #include <DrawTrSurf_Curve.hxx>
11 #include <DrawTrSurf_Curve2d.hxx>
12 #include <DrawTrSurf_BezierCurve.hxx>
13 #include <DrawTrSurf_BSplineCurve.hxx>
14 #include <DrawTrSurf_BezierCurve2d.hxx>
15 #include <DrawTrSurf_BSplineCurve2d.hxx>
16 #include <Draw_Marker3D.hxx>
17 #include <Draw_Marker2D.hxx>
18 #include <Draw.hxx>
19 #include <Draw_Interpretor.hxx>
20 #include <Draw_Color.hxx>
21 #include <Draw_Display.hxx>
22
23 #include <GeomAPI.hxx>
24 #include <GeomAPI_IntCS.hxx>
25 #include <GeomAPI_IntSS.hxx>
26
27 //#include <GeomLProp.hxx>
28 #include <GeomProjLib.hxx>
29 #include <BSplCLib.hxx>
30
31 #include <gp.hxx>
32 #include <gp_Pln.hxx>
33 #include <gp_Parab2d.hxx>
34 #include <gp_Elips2d.hxx>
35 #include <gp_Hypr2d.hxx>
36
37 #include <Geom_Line.hxx>
38 #include <Geom_Circle.hxx>
39 #include <Geom_Ellipse.hxx>
40 #include <Geom_Parabola.hxx>
41 #include <Geom_Hyperbola.hxx>
42 #include <Geom2d_Line.hxx>
43 #include <Geom2d_Circle.hxx>
44 #include <Geom2d_Ellipse.hxx>
45 #include <Geom2d_Parabola.hxx>
46 #include <Geom2d_Hyperbola.hxx>
47 #include <Geom2d_BSplineCurve.hxx>
48 #include <Geom2d_Curve.hxx>
49
50 #include <GccAna_Lin2dBisec.hxx>
51 #include <GccAna_Circ2dBisec.hxx>
52 #include <GccAna_CircLin2dBisec.hxx>
53 #include <GccAna_CircPnt2dBisec.hxx>
54 #include <GccAna_LinPnt2dBisec.hxx>
55 #include <GccAna_Pnt2dBisec.hxx>
56 #include <GccInt_Bisec.hxx>
57 #include <GccInt_IType.hxx>
58
59 #include <Geom_Plane.hxx>
60 #include <Geom_Curve.hxx>
61 #include <Geom2d_Curve.hxx>
62 #include <Geom2d_TrimmedCurve.hxx>
63 #include <Geom_TrimmedCurve.hxx>
64
65 #include <Law_BSpline.hxx>
66
67 #include <TColgp_Array1OfPnt.hxx>
68 #include <TColgp_Array1OfPnt2d.hxx>
69 #include <TColStd_Array1OfReal.hxx>
70 #include <TColStd_Array1OfInteger.hxx>
71
72 #include <Adaptor3d_HCurve.hxx>
73 #include <Adaptor3d_HSurface.hxx>
74 #include <Adaptor3d_CurveOnSurface.hxx>
75
76 #include <GeomAdaptor_HCurve.hxx>
77 #include <GeomAdaptor_HSurface.hxx>
78 #include <GeomAdaptor.hxx>
79 #include <Geom2dAdaptor_HCurve.hxx>
80
81 #include <GeomAbs_SurfaceType.hxx>
82 #include <GeomAbs_CurveType.hxx>
83
84 #include <ProjLib_CompProjectedCurve.hxx>
85 #include <ProjLib_HCompProjectedCurve.hxx>
86 #include <Approx_CurveOnSurface.hxx>
87 #include <Precision.hxx>
88 #include <Geom2dAdaptor.hxx>
89
90
91 #include <Precision.hxx>
92
93 #include <Geom_Surface.hxx>
94 #include <Adaptor2d_HCurve2d.hxx>
95 #include <stdio.h>
96 #include <BSplCLib.hxx>
97 #include <Geom_BSplineSurface.hxx>
98 #include <Geom_BSplineCurve.hxx>
99 #include <GCPnts_QuasiUniformDeflection.hxx>
100 #include <GCPnts_UniformDeflection.hxx>
101 #include <GCPnts_TangentialDeflection.hxx>
102 #include <GeomAPI_ExtremaCurveCurve.hxx>
103 #include <gce_MakeLin.hxx>
104 #include <TColStd_Array1OfBoolean.hxx>
105 #include <GeomAdaptor_HSurface.hxx>
106 #include <Adaptor3d_TopolTool.hxx>
107 #include <TColgp_Array2OfPnt.hxx>
108 #include <Geom_BSplineSurface.hxx>
109 #include <DrawTrSurf_BSplineSurface.hxx>
110 #include <TColStd_HArray1OfReal.hxx>
111
112 //epa test
113 #include <BRepBuilderAPI_MakeEdge.hxx>
114 #include <AIS_Shape.hxx>
115 #include <TopoDS_Edge.hxx>
116 #include <GeomLProp_CLProps.hxx>
117 #include <GCPnts_AbscissaPoint.hxx>
118 #include <GCPnts_UniformAbscissa.hxx>
119 #include <DBRep.hxx>
120
121 #ifdef WNT
122 Standard_IMPORT Draw_Viewer dout;
123 #endif
124
125 //=======================================================================
126 //function : polecurve2d
127 //purpose  : 
128 //=======================================================================
129
130 static Standard_Integer polelaw (Draw_Interpretor& , Standard_Integer n, const char** a)
131 {
132   Standard_Integer k,
133   jj,
134   qq,
135   i;
136
137
138   if (n < 3) return 1;
139   Standard_Boolean periodic = Standard_False ;
140   Standard_Integer deg = atoi(a[2]);
141   Standard_Integer nbk = atoi(a[3]);
142   
143   TColStd_Array1OfReal    knots(1, nbk);
144   TColStd_Array1OfInteger mults(1, nbk);
145   k = 4;
146   Standard_Integer Sigma = 0;
147   for (i = 1; i<=nbk; i++) {
148     knots( i) = atof(a[k]);
149     k++;
150     mults( i) = atoi(a[k]);
151     Sigma += mults(i);
152     k++;
153   }
154
155   Standard_Integer np;
156   np = Sigma - deg  -1;
157   TColStd_Array1OfReal flat_knots(1, Sigma) ;
158   jj = 1 ;
159   for (i = 1 ; i <= nbk ; i++) {
160     for(qq = 1 ; qq <= mults(i) ; qq++) {
161       flat_knots(jj) = knots(i) ;
162       jj ++ ;
163     }
164   }
165   
166   TColgp_Array1OfPnt2d poles  (1, np);
167   TColStd_Array1OfReal schoenberg_points(1,np) ;
168   BSplCLib::BuildSchoenbergPoints(deg,
169                                   flat_knots,
170                                   schoenberg_points) ;
171   for (i = 1; i <= np; i++) {
172     poles(i).SetCoord(schoenberg_points(i),atof(a[k]));
173     k++;
174   }
175     
176   Handle(Geom2d_BSplineCurve) result =
177     new Geom2d_BSplineCurve(poles, knots, mults, deg, periodic);
178   DrawTrSurf::Set(a[1],result);
179
180   
181   return 0;
182 }
183 //=======================================================================
184 //function : to2d
185 //purpose  : 
186 //=======================================================================
187
188 static Standard_Integer to2d (Draw_Interpretor& , Standard_Integer n, const char** a)
189 {
190   if (n < 3) return 1;
191
192   // get the curve
193   Handle(Geom_Curve) C = DrawTrSurf::GetCurve(a[2]);
194   if (C.IsNull())
195     return 1;
196
197   Handle(Geom_Surface) S;
198   if (n >= 4) {
199     S = DrawTrSurf::GetSurface(a[3]);
200     if (S.IsNull()) return 1;
201   }
202   else
203     S = new Geom_Plane(gp::XOY());
204   
205   Handle(Geom_Plane) P = Handle(Geom_Plane)::DownCast(S);
206   if (P.IsNull()) return 1;
207   Handle(Geom2d_Curve) r = GeomAPI::To2d(C,P->Pln());
208   DrawTrSurf::Set(a[1],r);
209   return 0;
210 }
211
212 //=======================================================================
213 //function : to3d
214 //purpose  : 
215 //=======================================================================
216
217 static Standard_Integer to3d (Draw_Interpretor& , Standard_Integer n, const char** a)
218 {
219   if (n < 3) return 1;
220   
221   Handle(Geom2d_Curve) C = DrawTrSurf::GetCurve2d(a[2]);
222   if (C.IsNull()) return 1;
223   
224   Handle(Geom_Surface) S;
225   if (n >= 4) {
226     S = DrawTrSurf::GetSurface(a[3]);
227     if (S.IsNull()) return 1;
228   }
229   else
230     S = new Geom_Plane(gp::XOY());
231   
232   Handle(Geom_Plane) P = Handle(Geom_Plane)::DownCast(S);
233   if (P.IsNull()) return 1;
234   Handle(Geom_Curve) r = GeomAPI::To3d(C,P->Pln());
235   
236   DrawTrSurf::Set(a[1],r);
237   return 0;
238 }
239
240 //=======================================================================
241 //function : gproject
242 //purpose  : 
243 //=======================================================================
244
245
246 static Standard_Integer gproject(Draw_Interpretor& di, Standard_Integer n, const char** a)
247 {
248   
249   char newname[1024];
250   char* temp = newname;
251   char newname1[10];
252   char* temp1 = newname1;
253   char name[100];
254   Standard_Integer ONE = 1;
255
256   if (n == 3)
257     sprintf(name,"p");
258   else if (n == 4) {
259     sprintf(name,"%s",a[1]);
260     ONE = 2;
261   }
262   else {
263    di << "gproject wait 2 or 3 arguments" << "\n";
264    return 1;
265   } 
266
267   Handle(Geom_Curve) Cur = DrawTrSurf::GetCurve(a[ONE]);
268   Handle(Geom_Surface) Sur = DrawTrSurf::GetSurface(a[ONE+1]);
269   if (Cur.IsNull() || Sur.IsNull()) return 1;
270
271   Handle(GeomAdaptor_HCurve) hcur = new GeomAdaptor_HCurve(Cur);
272   Handle(GeomAdaptor_HSurface) hsur = new GeomAdaptor_HSurface(Sur);
273
274
275   Standard_Real myTol3d = 1.e-6;
276   GeomAbs_Shape myContinuity = GeomAbs_C2;
277   Standard_Integer myMaxDegree = 14, myMaxSeg = 16;
278
279
280   ProjLib_CompProjectedCurve Projector(hsur, hcur, myTol3d/10, myTol3d/10);
281   Handle(ProjLib_HCompProjectedCurve) HProjector = new ProjLib_HCompProjectedCurve();
282   HProjector->Set(Projector);
283
284   Standard_Integer k;
285   Standard_Real Udeb, Ufin, UIso, VIso;
286   Standard_Integer Only2d, Only3d;
287   gp_Pnt2d P2d, Pdeb, Pfin;
288   gp_Pnt P;
289   Handle(Adaptor2d_HCurve2d) HPCur;
290   Handle(Geom2d_Curve) PCur2d; // Only for isoparametric projection
291
292   for(k = 1; k <= Projector.NbCurves(); k++){
293     sprintf(newname,"%s_%d",name,k);
294     sprintf(newname1,"%s2d_%d",name,k);
295     if(Projector.IsSinglePnt(k, P2d)){
296 //      cout<<"Part "<<k<<" of the projection is punctual"<<endl;
297       Projector.GetSurface()->D0(P2d.X(), P2d.Y(), P);
298       DrawTrSurf::Set(temp, P);
299       DrawTrSurf::Set(temp1, P2d);
300       di<<temp<<" is 3d projected curve"<<"\n";
301       di<<temp1<<" is pcurve"<<"\n";
302     }
303     else {
304       Only2d = Only3d = Standard_False;
305       Projector.Bounds(k, Udeb, Ufin);
306       gp_Dir2d Dir; // Only for isoparametric projection
307       
308       if (Projector.IsUIso(k, UIso)) {
309 //      cout<<"Part "<<k<<" of the projection is U-isoparametric curve"<<endl;
310         Projector.D0(Udeb, Pdeb);
311         Projector.D0(Ufin, Pfin);
312         Udeb = Pdeb.Y();
313         Ufin = Pfin.Y();
314         if (Udeb > Ufin) {
315           Dir = gp_Dir2d(0, -1);
316           Udeb = - Udeb;
317           Ufin = - Ufin;
318         }
319         else Dir = gp_Dir2d(0, 1);
320         PCur2d = new Geom2d_TrimmedCurve(new Geom2d_Line(gp_Pnt2d(UIso, 0), Dir), Udeb, Ufin);
321         HPCur = new Geom2dAdaptor_HCurve(PCur2d);
322         Only3d = Standard_True;
323       }
324       else if(Projector.IsVIso(k, VIso)) {
325 //      cout<<"Part "<<k<<" of the projection is V-isoparametric curve"<<endl;
326         Projector.D0(Udeb, Pdeb);
327         Projector.D0(Ufin, Pfin);
328         Udeb = Pdeb.X();
329         Ufin = Pfin.X();
330         if (Udeb > Ufin) {
331           Dir = gp_Dir2d(-1, 0);
332           Udeb = - Udeb;
333           Ufin = - Ufin;
334         }
335         else Dir = gp_Dir2d(1, 0);
336         PCur2d = new Geom2d_TrimmedCurve(new Geom2d_Line(gp_Pnt2d(0, VIso), Dir), Udeb, Ufin);
337         HPCur = new Geom2dAdaptor_HCurve(PCur2d);
338         Only3d = Standard_True;
339       }
340       else HPCur = HProjector;
341       
342       if(Projector.MaxDistance(k) <= myTol3d)
343         Only2d = Standard_True;
344       
345       if(Only2d && Only3d) {
346         Handle(Geom_Curve) OutCur = new Geom_TrimmedCurve(GeomAdaptor::MakeCurve(hcur->Curve()), Ufin, Udeb);
347         DrawTrSurf::Set(temp, OutCur);
348         DrawTrSurf::Set(temp1, PCur2d);
349         di<<temp<<" is 3d projected curve"<<"\n";
350         di<<temp1<<" is pcurve"<<"\n";
351         return 0;
352         }
353       else {
354         Approx_CurveOnSurface appr(HPCur, hsur, Udeb, Ufin, myTol3d, 
355                                    myContinuity, myMaxDegree, myMaxSeg, 
356                                    Only3d, Only2d);
357         if(!Only3d) {
358           PCur2d = appr.Curve2d();
359           di << " Error in 2d is " << appr.MaxError2dU()
360                << ";  "  << appr.MaxError2dV() << "\n"; 
361         }
362         if(Only2d) {
363           Handle(Geom_Curve) OutCur = 
364             new Geom_TrimmedCurve(GeomAdaptor::MakeCurve(hcur->Curve()), 
365                                   Ufin, Udeb);
366           DrawTrSurf::Set(temp, OutCur);
367           }
368         else {
369           di << " Error in 3d is " <<  appr.MaxError3d() << "\n";
370           DrawTrSurf::Set(temp, appr.Curve3d());
371         }
372         DrawTrSurf::Set(temp1, PCur2d);
373         di<<temp<<" is 3d projected curve"<<"\n";
374         di<<temp1<<" is pcurve"<<"\n";
375       }
376     }
377   }
378 return 0;
379 }
380 //=======================================================================
381 //function : project
382 //purpose  : 
383 //=======================================================================
384
385 static Standard_Integer project (Draw_Interpretor& di, 
386                                  Standard_Integer n, const char** a)
387 {
388   if ( n == 1) {
389     
390     di << "project result2d c3d surf [-e p] [-v n] [-t tol]" << "\n";
391     di << "   -e p   : extent the surface of <p>%" << "\n";
392     di << "   -v n   : verify the projection at <n> points." << "\n";
393     di << "   -t tol : set the tolerance for approximation" << "\n";
394     return 0;
395   }
396
397   if (n < 4) return 1;
398   Handle(Geom_Surface) GS = DrawTrSurf::GetSurface(a[3]);
399   if (GS.IsNull()) return 1;
400     
401   Handle(Geom_Curve) GC = DrawTrSurf::GetCurve(a[2]);
402   if (GC.IsNull()) return 1;
403
404   Standard_Real tolerance = Precision::Confusion() ;
405
406   Standard_Real U1,U2,V1,V2;
407   GS->Bounds(U1,U2,V1,V2);
408
409   Standard_Boolean Verif = Standard_False, Extent = Standard_False;
410   Standard_Integer NbPoints=0;
411
412   Standard_Integer index = 4;
413   while ( index+1 < n) {
414     if ( a[index][0] != '-') return 1;
415
416     if ( a[index][1] == 'e') {
417       Standard_Real p = atof(a[index+1]);
418       Standard_Real dU = p * (U2 - U1) / 100.;
419       Standard_Real dV = p * (V2 - V1) / 100.;
420       U1 -= dU; U2 += dU; V1 -= dV; V2 += dV;
421       Extent = Standard_True;
422     }
423     else if ( a[index][1] == 'v') {
424       Verif = Standard_True;
425       NbPoints = atoi(a[index+1]);
426     }
427     else if ( a[index][1] == 't') {
428       tolerance = atof(a[index+1]);
429     }
430     index += 2;
431   }
432   
433   Handle(Geom2d_Curve) G2d = 
434     GeomProjLib::Curve2d(GC, GS, U1, U2, V1, V2, tolerance);
435
436   if ( G2d.IsNull() ) {
437     di << "\n" << "Projection Failed" << "\n";
438     return 1;
439   }
440   else {
441     DrawTrSurf::Set(a[1],G2d);
442   }
443   if ( Verif) {  // verify the projection on n points
444     if ( NbPoints <= 0) { 
445       di << " n must be positive" << "\n";
446       return 0;
447     }
448     gp_Pnt P1,P2;
449     gp_Pnt2d P2d;
450     
451     Standard_Real U, dU;
452     Standard_Real Dist,DistMax = -1.;
453     U1 = GC->FirstParameter();
454     U2 = GC->LastParameter();
455     dU = ( U2 - U1) / (NbPoints + 1);
456     for ( Standard_Integer i = 0 ; i <= NbPoints +1; i++) {
457       U = U1 + i *dU;
458       P1 = GC->Value(U);
459       P2d = G2d->Value(U);
460       P2 = GS->Value(P2d.X(), P2d.Y());
461       Dist = P1.Distance(P2);
462       di << " Parameter = " << U << "\tDistance = " << Dist << "\n";
463       if ( Dist > DistMax) DistMax = Dist;
464     }
465     di << " **** Distance Maximale : " << DistMax << "\n";
466   }
467
468   return 0;
469 }
470
471 //=======================================================================
472 //function : projonplane
473 //purpose  : 
474 //=======================================================================
475
476 Standard_Integer projonplane(Draw_Interpretor& di, 
477                              Standard_Integer n, const char** a)
478 {
479   if ( n < 4 ) return 1;
480
481   Handle(Geom_Surface) S = DrawTrSurf::GetSurface(a[3]);
482   if ( S.IsNull()) return 1;
483
484   Handle(Geom_Plane)   Pl = Handle(Geom_Plane)::DownCast(S);
485   if ( Pl.IsNull()) {
486     di << " The surface must be a plane" << "\n";
487     return 1;
488   }
489
490   Handle(Geom_Curve) C = DrawTrSurf::GetCurve(a[2]);
491   if ( C.IsNull()) return 1;
492   
493   Standard_Boolean Param = Standard_True;
494   if ((n == 5 && atoi(a[4]) == 0) ||
495       (n == 8 && atoi(a[7]) == 0)) Param = Standard_False;
496
497   gp_Dir D;
498   
499   if ( n == 8) {
500     D = gp_Dir(atof(a[4]),atof(a[5]),atof(a[6]));
501   }
502   else { 
503     D = Pl->Pln().Position().Direction();
504   }
505
506   Handle(Geom_Curve) Res = 
507     GeomProjLib::ProjectOnPlane(C,Pl,D,Param);
508
509   DrawTrSurf::Set(a[1],Res);
510   return 0;
511
512 }
513
514
515 //=======================================================================
516 //function : bisec
517 //purpose  : 
518 //=======================================================================
519
520 static void solution(const Handle(GccInt_Bisec)& Bis,
521                      const char* name,
522                      const Standard_Integer i)
523 {
524   char solname[200];
525   if ( i == 0) 
526     sprintf(solname,"%s",name);
527   else
528     sprintf(solname,"%s_%d",name,i);
529   const char* temp = solname; // pour portage WNT
530
531   switch ( Bis->ArcType()) {
532   case GccInt_Lin:
533     DrawTrSurf::Set(temp, new Geom2d_Line(Bis->Line()));
534     break;
535   case GccInt_Cir:
536     DrawTrSurf::Set(temp, new Geom2d_Circle(Bis->Circle()));
537     break;
538   case GccInt_Ell:
539     DrawTrSurf::Set(temp, new Geom2d_Ellipse(Bis->Ellipse()));
540     break;
541   case GccInt_Par:
542     DrawTrSurf::Set(temp, new Geom2d_Parabola(Bis->Parabola()));
543     break;
544   case GccInt_Hpr:
545     DrawTrSurf::Set(temp, new Geom2d_Hyperbola(Bis->Hyperbola()));
546     break;
547   case GccInt_Pnt:
548     DrawTrSurf::Set(temp, Bis->Point());
549     break;
550   }
551 }
552
553 static Standard_Integer bisec (Draw_Interpretor& di, 
554                                Standard_Integer n, const char** a)
555 {
556   if (n < 4) return 1;
557   
558   Handle(Geom2d_Curve) C1 = DrawTrSurf::GetCurve2d(a[2]);
559   Handle(Geom2d_Curve) C2 = DrawTrSurf::GetCurve2d(a[3]);
560   gp_Pnt2d P1,P2;
561   Standard_Boolean ip1 = DrawTrSurf::GetPoint2d(a[2],P1);
562   Standard_Boolean ip2 = DrawTrSurf::GetPoint2d(a[3],P2);
563   Standard_Integer i, Compt = 0;
564   Standard_Integer NbSol = 0;
565
566   if ( !C1.IsNull()) {
567     Handle(Standard_Type) Type1 = C1->DynamicType();
568     if ( !C2.IsNull()) {
569       Handle(Standard_Type) Type2 = C2->DynamicType();
570       if ( Type1 == STANDARD_TYPE(Geom2d_Line) && 
571            Type2 == STANDARD_TYPE(Geom2d_Line)   ) {
572         GccAna_Lin2dBisec Bis(Handle(Geom2d_Line)::DownCast(C1)->Lin2d(),
573                               Handle(Geom2d_Line)::DownCast(C2)->Lin2d());
574         if ( Bis.IsDone()) {
575           char solname[200];
576           NbSol = Bis.NbSolutions();
577           for ( i = 1; i <= NbSol; i++) {
578             sprintf(solname,"%s_%d",a[1],i);
579             const char* temp = solname; // pour portage WNT
580             DrawTrSurf::Set(temp,new Geom2d_Line(Bis.ThisSolution(i)));
581           }
582         }
583         else {
584           di << " Bisec has failed !!" << "\n";
585           return 1;
586         }
587       }
588       else if ( Type1 == STANDARD_TYPE(Geom2d_Line) && 
589                 Type2 == STANDARD_TYPE(Geom2d_Circle) ) {
590         GccAna_CircLin2dBisec 
591           Bis(Handle(Geom2d_Circle)::DownCast(C2)->Circ2d(),
592               Handle(Geom2d_Line)::DownCast(C1)->Lin2d());
593         if ( Bis.IsDone()) {
594           NbSol= Bis.NbSolutions();
595           if ( NbSol >= 2) Compt = 1;
596           for ( i = 1; i <= NbSol; i++) {
597             solution(Bis.ThisSolution(i),a[1],Compt);
598             Compt++;
599           }
600         }
601         else {
602           di << " Bisec has failed !!" << "\n";
603           return 1;
604         }
605       }
606       else if ( Type2 == STANDARD_TYPE(Geom2d_Line) && 
607                 Type1 == STANDARD_TYPE(Geom2d_Circle) ) {
608         GccAna_CircLin2dBisec 
609           Bis(Handle(Geom2d_Circle)::DownCast(C1)->Circ2d(), 
610               Handle(Geom2d_Line)::DownCast(C2)->Lin2d());
611         if ( Bis.IsDone()) {
612 //        char solname[200];
613           NbSol = Bis.NbSolutions();
614           if ( NbSol >= 2) Compt = 1;
615           for ( i = 1; i <= NbSol; i++) {
616             solution(Bis.ThisSolution(i),a[1],Compt);
617             Compt++;
618           }
619         }
620         else {
621           di << " Bisec has failed !!" << "\n";
622           return 1;
623         }
624       }
625       else if ( Type2 == STANDARD_TYPE(Geom2d_Circle) && 
626                 Type1 == STANDARD_TYPE(Geom2d_Circle) ) {
627         GccAna_Circ2dBisec 
628           Bis(Handle(Geom2d_Circle)::DownCast(C1)->Circ2d(), 
629               Handle(Geom2d_Circle)::DownCast(C2)->Circ2d());
630         if ( Bis.IsDone()) {
631 //        char solname[200];
632           NbSol = Bis.NbSolutions();
633           if ( NbSol >= 2) Compt = 1;
634           for ( i = 1; i <= NbSol; i++) {
635             solution(Bis.ThisSolution(i),a[1],Compt);
636             Compt++;
637           }
638         }
639         else {
640           di << " Bisec has failed !!" << "\n";
641           return 1;
642         }
643       }
644       else {
645         di << " args must be line/circle/point line/circle/point" << "\n";
646         return 1;
647       }
648     }
649     else if (ip2) {
650       if ( Type1 == STANDARD_TYPE(Geom2d_Circle)) {
651         GccAna_CircPnt2dBisec Bis
652           (Handle(Geom2d_Circle)::DownCast(C1)->Circ2d(),P2);
653         if ( Bis.IsDone()) {
654           NbSol = Bis.NbSolutions();
655           if ( NbSol >= 2) Compt = 1;
656           for ( i = 1; i <= NbSol; i++) {
657             solution(Bis.ThisSolution(i),a[1],Compt);
658             Compt++;
659           }
660         }
661         else {
662           di << " Bisec has failed !!" << "\n";
663           return 1;
664         }
665       }
666       else if ( Type1 == STANDARD_TYPE(Geom2d_Line)) {
667         GccAna_LinPnt2dBisec Bis
668           (Handle(Geom2d_Line)::DownCast(C1)->Lin2d(),P2);
669         if ( Bis.IsDone()) {
670           NbSol = 1;
671           solution(Bis.ThisSolution(),a[1],0);
672         }
673         else {
674           di << " Bisec has failed !!" << "\n";
675           return 1;
676         }
677       }
678     }
679     else {
680       di << " the second arg must be line/circle/point " << "\n";
681     }
682   }
683   else if ( ip1) {
684     if ( !C2.IsNull()) {
685       Handle(Standard_Type) Type2 = C2->DynamicType();
686       if ( Type2 == STANDARD_TYPE(Geom2d_Circle)) {
687         GccAna_CircPnt2dBisec Bis
688           (Handle(Geom2d_Circle)::DownCast(C2)->Circ2d(),P1);
689         if ( Bis.IsDone()) {
690           NbSol = Bis.NbSolutions();
691           if ( NbSol >= 2) Compt = 1;
692           for ( i = 1; i <= Bis.NbSolutions(); i++) {
693             solution(Bis.ThisSolution(i),a[1],Compt);
694             Compt++;
695           }
696         }
697         else {
698           di << " Bisec has failed !!" << "\n";
699           return 1;
700         }
701       }
702       else if ( Type2 == STANDARD_TYPE(Geom2d_Line)) {
703         GccAna_LinPnt2dBisec Bis
704           (Handle(Geom2d_Line)::DownCast(C2)->Lin2d(),P1);
705         if ( Bis.IsDone()) {
706           NbSol = 1;
707           solution(Bis.ThisSolution(),a[1],0);
708         }
709         else {
710           di << " Bisec has failed !!" << "\n";
711           return 1;
712         }
713       }
714     }
715     else if (ip2) {
716       GccAna_Pnt2dBisec Bis(P1,P2);
717       if ( Bis.HasSolution()) {
718         NbSol = 1;
719         DrawTrSurf::Set(a[1],new Geom2d_Line(Bis.ThisSolution()));
720       }
721       else {
722         di << " Bisec has failed !!" << "\n";
723         return 1;
724       }
725     }
726     else {
727       di << " the second arg must be line/circle/point " << "\n";
728       return 1;
729     }
730   }
731   else {
732     di << " args must be line/circle/point line/circle/point" << "\n";
733     return 1;
734   }
735
736   if ( NbSol >= 2) {
737     di << "There are " << NbSol << " Solutions." << "\n";
738   }
739   else {
740     di << "There is " << NbSol << " Solution." << "\n";
741   }
742
743   return 0;
744 }
745
746 //=======================================================================
747 //function : cmovetangent
748 //purpose  : 
749 //=======================================================================
750
751 static Standard_Integer movelaw (Draw_Interpretor& di, Standard_Integer n, const char** a)
752 {
753   Standard_Integer dimension,
754   ii,
755   condition=0,
756   error_status ;
757   Standard_Real u,
758   x,
759   tolerance,
760   tx ;
761
762   u = atof(a[2]);
763   x = atof(a[3]);
764   tolerance = 1.0e-5 ;
765   dimension = 2 ;
766   if (n < 5) {
767       return 1 ;
768   }
769   Handle(Geom2d_BSplineCurve) G2 = DrawTrSurf::GetBSplineCurve2d(a[1]);
770   if (!G2.IsNull()) {
771       tx = atof(a[4]) ;
772       if (n == 6) {
773         condition = Max(atoi(a[5]), -1)  ;
774         condition = Min(condition, G2->Degree()-1) ;
775       }
776       TColgp_Array1OfPnt2d   curve_poles(1,G2->NbPoles()) ;
777       TColStd_Array1OfReal    law_poles(1,G2->NbPoles()) ;
778       TColStd_Array1OfReal    law_knots(1,G2->NbKnots()) ;
779       TColStd_Array1OfInteger law_mults(1,G2->NbKnots()) ;
780
781       G2->Knots(law_knots) ;
782       G2->Multiplicities(law_mults) ;
783       G2->Poles(curve_poles) ;
784       for (ii = 1 ; ii <= G2->NbPoles() ; ii++) {
785         law_poles(ii) = curve_poles(ii).Coord(2) ;
786       }
787
788       Law_BSpline  a_law(law_poles,
789                          law_knots,
790                          law_mults,
791                          G2->Degree(),
792                          Standard_False) ;
793
794       a_law.MovePointAndTangent(u, 
795                               x, 
796                               tx,
797                               tolerance,
798                               condition,
799                               condition,
800                               error_status) ;
801
802       for (ii = 1 ; ii <= G2->NbPoles() ; ii++) {
803         curve_poles(ii).SetCoord(2,a_law.Pole(ii)) ;
804         G2->SetPole(ii,curve_poles(ii)) ;
805       }
806
807         
808       if (! error_status) {
809         Draw::Repaint();
810         }
811       else {
812         di << "Not enought degree of freedom increase degree please" << "\n";
813       }
814       
815     
816     }
817   return 0;
818
819
820
821 //Static method computing deviation of curve and polyline
822
823 static void ComputeDeviation(const Handle(Geom_Curve)& theCurve,
824                              const Handle(Geom_BSplineCurve)& thePnts,
825                              Standard_Real& theDmax,
826                              Standard_Real& theUfMax,
827                              Standard_Real& theUlMax,
828                              Standard_Integer& theImax)
829 {
830   theDmax = 0.;
831   theUfMax = 0.;
832   theUlMax = 0.;
833   theImax = 0;
834   Standard_Real dmax = 0., ufmax = 0., ulmax = 0.;
835   Standard_Integer imax = 0;
836
837   //take knots
838   Standard_Integer nbp = thePnts->NbKnots();
839   TColStd_Array1OfReal aKnots(1, nbp);
840   thePnts->Knots(aKnots);
841
842   Standard_Integer i;
843   for(i = 1; i < nbp; ++i) {
844     Standard_Real uf = aKnots(i);
845     Standard_Real ul = aKnots(i+1);
846
847     GeomAPI_ExtremaCurveCurve ECC(theCurve, thePnts, uf, ul, uf, ul);
848
849     Standard_Integer nbe = ECC.NbExtrema();
850     if(nbe > 0) {
851       Standard_Integer k;
852       Standard_Real d = 0.;
853       for(k = 1; k <= nbe; k++) {
854         if(ECC.Distance(k) > d) d = ECC.Distance(k);
855       }
856
857       if(d > theDmax) {
858         theDmax = d;
859               theUfMax = uf;
860               theUlMax = ul;
861               theImax = i;
862       }
863     }
864   }
865 }
866
867
868 //=======================================================================
869 //function : crvpoints
870 //purpose  : 
871 //=======================================================================
872
873 static Standard_Integer crvpoints (Draw_Interpretor& di, Standard_Integer /*n*/, const char** a)
874 {
875   Standard_Integer i, nbp;
876   Standard_Real defl;
877
878   Handle(Geom_Curve) C = DrawTrSurf::GetCurve(a[2]);
879   defl = atof(a[3]);
880
881   GeomAdaptor_Curve GAC(C);
882   GCPnts_QuasiUniformDeflection PntGen(GAC, defl);
883     
884   if(!PntGen.IsDone()) {
885     di << "Points generation failed" << "\n";
886     return 1;
887   }
888
889   nbp = PntGen.NbPoints();
890   di << "Nb points : " << nbp << "\n";
891
892   TColgp_Array1OfPnt aPoles(1, nbp);
893   TColStd_Array1OfReal aKnots(1, nbp);
894   TColStd_Array1OfInteger aMults(1, nbp);
895
896   for(i = 1; i <= nbp; ++i) {
897     aPoles(i) = PntGen.Value(i);
898     aKnots(i) = PntGen.Parameter(i);
899     aMults(i) = 1;
900   }
901   
902   aMults(1) = 2;
903   aMults(nbp) = 2;
904
905   Handle(Geom_BSplineCurve) aPnts = new Geom_BSplineCurve(aPoles, aKnots, aMults, 1);
906   Handle(DrawTrSurf_BSplineCurve) aDrCrv = new DrawTrSurf_BSplineCurve(aPnts);
907
908   aDrCrv->ClearPoles();
909   Draw_Color aKnColor(Draw_or);
910   aDrCrv->SetKnotsColor(aKnColor);
911   aDrCrv->SetKnotsShape(Draw_Plus);
912
913   Draw::Set(a[1], aDrCrv);
914
915   Standard_Real dmax = 0., ufmax = 0., ulmax = 0.;
916   Standard_Integer imax = 0;
917
918   //check deviation
919   ComputeDeviation(C,aPnts,dmax,ufmax,ulmax,imax);
920   di << "Max defl: " << dmax << " " << ufmax << " " << ulmax << " " << i << "\n"; 
921
922   return 0;
923
924
925 //=======================================================================
926 //function : crvtpoints
927 //purpose  : 
928 //=======================================================================
929
930 static Standard_Integer crvtpoints (Draw_Interpretor& di, Standard_Integer n, const char** a)
931 {
932   Standard_Integer i, nbp;
933   Standard_Real defl, angle = Precision::Angular();
934
935   Handle(Geom_Curve) C = DrawTrSurf::GetCurve(a[2]);
936   defl = atof(a[3]);
937
938   if(n > 3)
939     angle = atof(a[4]);
940
941   GeomAdaptor_Curve GAC(C);
942   GCPnts_TangentialDeflection PntGen(GAC, angle, defl, 2);
943   
944   nbp = PntGen.NbPoints();
945   di << "Nb points : " << nbp << "\n";
946
947   TColgp_Array1OfPnt aPoles(1, nbp);
948   TColStd_Array1OfReal aKnots(1, nbp);
949   TColStd_Array1OfInteger aMults(1, nbp);
950
951   for(i = 1; i <= nbp; ++i) {
952     aPoles(i) = PntGen.Value(i);
953     aKnots(i) = PntGen.Parameter(i);
954     aMults(i) = 1;
955   }
956   
957   aMults(1) = 2;
958   aMults(nbp) = 2;
959
960   Handle(Geom_BSplineCurve) aPnts = new Geom_BSplineCurve(aPoles, aKnots, aMults, 1);
961   Handle(DrawTrSurf_BSplineCurve) aDrCrv = new DrawTrSurf_BSplineCurve(aPnts);
962
963   aDrCrv->ClearPoles();
964   Draw_Color aKnColor(Draw_or);
965   aDrCrv->SetKnotsColor(aKnColor);
966   aDrCrv->SetKnotsShape(Draw_Plus);
967
968   Draw::Set(a[1], aDrCrv);
969
970   Standard_Real dmax = 0., ufmax = 0., ulmax = 0.;
971   Standard_Integer imax = 0;
972
973   //check deviation
974   ComputeDeviation(C,aPnts,dmax,ufmax,ulmax,imax);
975   di << "Max defl: " << dmax << " " << ufmax << " " << ulmax << " " << i << "\n"; 
976
977   return 0;
978
979 //=======================================================================
980 //function : uniformAbscissa
981 //purpose  : epa test (TATA-06-002 (Problem with GCPnts_UniformAbscissa class)
982 //=======================================================================
983 static Standard_Integer uniformAbscissa (Draw_Interpretor& di, Standard_Integer n, const char** a)
984 {
985   if( n != 3 )
986     return 1;
987   
988   /*Handle(Geom_BSplineCurve) ellip;
989   ellip = DrawTrSurf::GetBSplineCurve(a[1]);
990   if (ellip.IsNull())
991   {
992     di << " BSpline is NULL  "<<"\n";     
993     return 1;
994   }*/
995
996   Handle(Geom_Curve) ellip;
997   ellip = DrawTrSurf::GetCurve(a[1]);
998   if (ellip.IsNull())
999   {
1000     di << " Curve is NULL  "<<"\n";     
1001     return 1;
1002   }
1003
1004   Standard_Integer nocp;
1005   nocp = atoi(a[2]);
1006   if(nocp < 2)
1007     return 1;
1008
1009
1010   //test nbPoints for Geom_Ellipse
1011
1012   try
1013   {
1014     GeomLProp_CLProps Prop(ellip,2,Precision::Intersection());
1015     Prop.SetCurve(ellip);
1016
1017     GeomAdaptor_Curve GAC(ellip);
1018     di<<"Type Of curve: "<<GAC.GetType()<<"\n";
1019     Standard_Real Tol = Precision::Confusion();
1020     Standard_Real L;
1021
1022     L = GCPnts_AbscissaPoint::Length(GAC, GAC.FirstParameter(), GAC.LastParameter(), Tol);
1023     di<<"Ellipse length = "<<L<<"\n";
1024     Standard_Real Abscissa = L/(nocp-1);
1025     di << " CUR : Abscissa " << Abscissa << "\n";
1026
1027     GCPnts_UniformAbscissa myAlgo(GAC, Abscissa, ellip->FirstParameter(), ellip->LastParameter());
1028     if ( myAlgo.IsDone() )
1029     {
1030       di << " CasCurve  - nbpoints " << myAlgo.NbPoints() << "\n";
1031       for(Standard_Integer i = 1; i<= myAlgo.NbPoints(); i++ )
1032         di << i <<" points = " << myAlgo.Parameter( i ) << "\n";
1033     }
1034   }
1035
1036   catch (Standard_Failure )
1037   {
1038     di << " Standard Failure  " <<"\n";
1039   }
1040   return 0;
1041 }
1042
1043 //=======================================================================
1044 //function : EllipsUniformAbscissa
1045 //purpose  : epa test (TATA-06-002 (Problem with GCPnts_UniformAbscissa class)
1046 //=======================================================================
1047 static Standard_Integer EllipsUniformAbscissa (Draw_Interpretor& di, Standard_Integer n, const char** a)
1048 {
1049   if( n != 4 )
1050     return 1;  
1051   
1052   Standard_Real R1;
1053   R1 = atof(a[1]);
1054   Standard_Real R2;
1055   R2 = atof(a[2]);
1056
1057   Standard_Integer nocp;
1058   nocp = atoi(a[3]);
1059   if(nocp < 2)
1060     return 1;
1061   
1062   //test nbPoints for Geom_Ellipse
1063   Handle_Geom_Ellipse ellip;
1064
1065
1066   try
1067   {
1068     gp_Pnt location;
1069     location = gp_Pnt( 0.0, 0.0, 0.0);
1070     gp_Dir main_direction(0.0, 0.0, 1.0);
1071
1072     gp_Dir x_direction(1.0, 0.0, 0.0);
1073     gp_Ax2 mainaxis( location, main_direction);
1074
1075     mainaxis.SetXDirection(x_direction);
1076     ellip = new Geom_Ellipse(mainaxis,R1, R2);
1077
1078     BRepBuilderAPI_MakeEdge curve_edge(ellip);
1079     TopoDS_Edge edge_curve = curve_edge.Edge();
1080
1081     DBRep::Set("Ellipse",edge_curve);
1082   }
1083   
1084   catch(Standard_Failure)
1085   {
1086     di << " Standard Failure  "<<"\n";     
1087   }
1088
1089   try
1090   {
1091     GeomLProp_CLProps Prop(ellip,2,Precision::Intersection());
1092     Prop.SetCurve(ellip);
1093
1094     GeomAdaptor_Curve GAC(ellip);
1095     di<<"Type Of curve: "<<GAC.GetType()<<"\n";
1096     Standard_Real Tol = Precision::Confusion();
1097     Standard_Real L;
1098
1099     L = GCPnts_AbscissaPoint::Length(GAC, GAC.FirstParameter(), GAC.LastParameter(), Tol);
1100     di<<"Ellipse length = "<<L<<"\n";
1101     Standard_Real Abscissa = L/(nocp-1);
1102     di << " CUR : Abscissa " << Abscissa << "\n";
1103
1104     GCPnts_UniformAbscissa myAlgo(GAC, Abscissa, ellip->FirstParameter(), ellip->LastParameter());
1105     if ( myAlgo.IsDone() )
1106     {
1107       di << " CasCurve  - nbpoints " << myAlgo.NbPoints() << "\n";
1108       for(Standard_Integer i = 1; i<= myAlgo.NbPoints(); i++ )
1109         di << i <<" points = " << myAlgo.Parameter( i ) << "\n";
1110     }
1111   }
1112
1113   catch (Standard_Failure )
1114   {
1115     di << " Standard Failure  " <<"\n";
1116   }
1117   return 0;
1118 }
1119
1120 //=======================================================================
1121 //function : mypoints
1122 //purpose  : 
1123 //=======================================================================
1124
1125 static Standard_Integer mypoints (Draw_Interpretor& di, Standard_Integer /*n*/, const char** a)
1126 {
1127   Standard_Integer i, nbp;
1128   Standard_Real defl;
1129
1130   Handle(Geom_Curve) C = DrawTrSurf::GetCurve(a[2]);
1131   defl = atof(a[3]);
1132   const Handle(Geom_BSplineCurve)& aBS = Handle(Geom_BSplineCurve)::DownCast(C);
1133
1134   if(aBS.IsNull()) return 1;
1135
1136   Standard_Integer ui1 = aBS->FirstUKnotIndex();
1137   Standard_Integer ui2 = aBS->LastUKnotIndex();
1138
1139   Standard_Integer nbsu = ui2-ui1+1; nbsu += (nbsu - 1) * (aBS->Degree()-1);
1140
1141   TColStd_Array1OfReal anUPars(1, nbsu);
1142   TColStd_Array1OfBoolean anUFlg(1, nbsu);
1143
1144   Standard_Integer j, k, nbi;
1145   Standard_Real t1, t2, dt;
1146
1147   //Filling of sample parameters
1148   nbi = aBS->Degree();
1149   k = 0;
1150   t1 = aBS->Knot(ui1);
1151   for(i = ui1+1; i <= ui2; ++i) {
1152     t2 = aBS->Knot(i);
1153     dt = (t2 - t1)/nbi;
1154     j = 1;
1155     do { 
1156       ++k;
1157       anUPars(k) = t1;
1158       anUFlg(k) = Standard_False;
1159       t1 += dt; 
1160     }
1161     while (++j <= nbi);
1162     t1 = t2;
1163   }
1164   ++k;
1165   anUPars(k) = t1;
1166
1167   Standard_Integer l;
1168   defl *= defl;
1169
1170   j = 1;
1171   anUFlg(1) = Standard_True;
1172   anUFlg(nbsu) = Standard_True;
1173   Standard_Boolean bCont = Standard_True;
1174   while (j < nbsu-1 && bCont) {
1175     t2 = anUPars(j);
1176     gp_Pnt p1 = aBS->Value(t2);
1177     for(k = j+2; k <= nbsu; ++k) {
1178       t2 = anUPars(k);
1179       gp_Pnt p2 = aBS->Value(t2);
1180       gce_MakeLin MkLin(p1, p2);
1181       const gp_Lin& lin = MkLin.Value();
1182       Standard_Boolean ok = Standard_True;
1183       for(l = j+1; l < k; ++l) {
1184         if(anUFlg(l)) continue;
1185         gp_Pnt pp =  aBS->Value(anUPars(l));
1186         Standard_Real d = lin.SquareDistance(pp);
1187           
1188         if(d <= defl) continue;
1189
1190         ok = Standard_False;
1191         break;
1192       }
1193
1194
1195       if(!ok) {
1196         j = k - 1;
1197         anUFlg(j) = Standard_True;
1198         break;
1199       }
1200
1201     }
1202
1203     if(k >= nbsu) bCont = Standard_False;
1204   }
1205
1206   nbp = 0;
1207   for(i = 1; i <= nbsu; ++i) {
1208     if(anUFlg(i)) nbp++;
1209   }
1210
1211   TColgp_Array1OfPnt aPoles(1, nbp);
1212   TColStd_Array1OfReal aKnots(1, nbp);
1213   TColStd_Array1OfInteger aMults(1, nbp);
1214   j = 0;
1215   for(i = 1; i <= nbsu; ++i) {
1216     if(anUFlg(i)) {
1217       ++j;
1218       aKnots(j) = anUPars(i);
1219       aMults(j) = 1;
1220       aPoles(j) = aBS->Value(aKnots(j));
1221     }
1222   }
1223   
1224   aMults(1) = 2;
1225   aMults(nbp) = 2;
1226
1227   Handle(Geom_BSplineCurve) aPnts = new Geom_BSplineCurve(aPoles, aKnots, aMults, 1);
1228   Handle(DrawTrSurf_BSplineCurve) aDrCrv = new DrawTrSurf_BSplineCurve(aPnts);
1229
1230   aDrCrv->ClearPoles();
1231   Draw_Color aKnColor(Draw_or);
1232   aDrCrv->SetKnotsColor(aKnColor);
1233   aDrCrv->SetKnotsShape(Draw_Plus);
1234
1235   Draw::Set(a[1], aDrCrv);
1236
1237   Standard_Real dmax = 0., ufmax = 0., ulmax = 0., uf, ul;
1238   Standard_Integer imax = 0;
1239
1240   ComputeDeviation(C,aPnts,dmax,ufmax,ulmax,imax);
1241   di << "Max defl: " << dmax << " " << ufmax << " " << ulmax << " " << imax << "\n"; 
1242
1243   return 0;
1244
1245
1246
1247
1248 //=======================================================================
1249 //function : surfpoints
1250 //purpose  : 
1251 //=======================================================================
1252
1253 static Standard_Integer surfpoints (Draw_Interpretor& /*di*/, Standard_Integer /*n*/, const char** a)
1254 {
1255   Standard_Integer i;
1256   Standard_Real defl;
1257
1258   Handle(Geom_Surface) S = DrawTrSurf::GetSurface(a[2]);
1259   defl = atof(a[3]);
1260
1261   Handle(GeomAdaptor_HSurface) AS = new GeomAdaptor_HSurface(S);
1262
1263   Handle(Adaptor3d_TopolTool) aTopTool = new Adaptor3d_TopolTool(AS);
1264
1265   aTopTool->SamplePnts(defl, 10, 10);
1266
1267   Standard_Integer nbpu = aTopTool->NbSamplesU();
1268   Standard_Integer nbpv = aTopTool->NbSamplesV();
1269   TColStd_Array1OfReal Upars(1, nbpu), Vpars(1, nbpv);
1270   aTopTool->UParameters(Upars);
1271   aTopTool->VParameters(Vpars);
1272
1273   TColgp_Array2OfPnt aPoles(1, nbpu, 1, nbpv);
1274   TColStd_Array1OfReal anUKnots(1, nbpu);
1275   TColStd_Array1OfReal aVKnots(1, nbpv);
1276   TColStd_Array1OfInteger anUMults(1, nbpu);
1277   TColStd_Array1OfInteger aVMults(1, nbpv);
1278
1279   Standard_Integer j;
1280   for(i = 1; i <= nbpu; ++i) {
1281     anUKnots(i) = Upars(i);
1282     anUMults(i) = 1;
1283     for(j = 1; j <= nbpv; ++j) {
1284       aVKnots(j) = Vpars(j);
1285       aVMults(j) = 1;
1286       aPoles(i,j) = S->Value(anUKnots(i),aVKnots(j));
1287     }
1288   }
1289   
1290   anUMults(1) = 2;
1291   anUMults(nbpu) = 2;
1292   aVMults(1) = 2;
1293   aVMults(nbpv) = 2;
1294
1295   Handle(Geom_BSplineSurface) aPnts = new Geom_BSplineSurface(aPoles, anUKnots,  aVKnots, 
1296                                                               anUMults, aVMults, 1, 1);
1297   Handle(DrawTrSurf_BSplineSurface) aDrSurf = new DrawTrSurf_BSplineSurface(aPnts);
1298
1299   aDrSurf->ClearPoles();
1300   Draw_Color aKnColor(Draw_or);
1301   aDrSurf->SetKnotsColor(aKnColor);
1302   aDrSurf->SetKnotsShape(Draw_Plus);
1303
1304   Draw::Set(a[1], aDrSurf);
1305
1306
1307   return 0;
1308
1309
1310
1311
1312 //=======================================================================
1313 //function : intersect
1314 //purpose  : 
1315 //=======================================================================
1316 static Standard_Integer intersection (Draw_Interpretor& di, Standard_Integer n, const char** a)
1317 {
1318   if (n < 4) {
1319     return 1;
1320   }
1321   //
1322   Handle(Geom_Curve) GC1;
1323   Handle(Geom_Surface) GS1 = DrawTrSurf::GetSurface(a[2]);
1324   if (GS1.IsNull()) {
1325     GC1 = DrawTrSurf::GetCurve(a[2]);
1326     if (GC1.IsNull())
1327       return 1;
1328   }
1329   //
1330   Handle(Geom_Surface) GS2 = DrawTrSurf::GetSurface(a[3]);
1331   if (GS2.IsNull()) {
1332     return 1;
1333   }
1334   //
1335   Standard_Real tol = Precision::Confusion();
1336   if (n == 5 || n == 9 || n == 13 || n == 17) tol = atof(a[n-1]);
1337   //
1338   Handle(Geom_Curve) Result;
1339   gp_Pnt             Point;
1340   //
1341   if (GC1.IsNull()) {
1342     GeomInt_IntSS Inters;
1343     //
1344     // Surface Surface
1345     if (n <= 5) {
1346       // General case
1347       Inters.Perform(GS1,GS2,tol,Standard_True);
1348     }
1349     else if (n == 8 || n == 9 || n == 12 || n == 13 || n == 16 || n == 17) {
1350       Standard_Boolean useStart = Standard_True, useBnd = Standard_True;
1351       Standard_Integer ista1=0,ista2=0,ibnd1=0,ibnd2=0;
1352       Standard_Real UVsta[4];
1353       Handle(GeomAdaptor_HSurface) AS1,AS2;
1354       //
1355       if (n <= 9) {        // user starting point
1356         useBnd = Standard_False;
1357         ista1 = 4; ista2 = 7;
1358       }
1359       else if (n <= 13) {        // user bounding
1360         useStart = Standard_False;
1361         ibnd1 = 4; ibnd2 = 11;
1362       }
1363       else {        // both user starting point and bounding
1364         ista1 = 4; ista2 = 7;
1365         ibnd1 = 8; ibnd2 = 15;
1366       }
1367       if (useStart)
1368         for (Standard_Integer i=ista1; i <= ista2; i++)
1369           UVsta[i-ista1] = atof(a[i]);
1370       if (useBnd) {
1371         Standard_Real UVbnd[8];
1372         for (Standard_Integer i=ibnd1; i <= ibnd2; i++)
1373           UVbnd[i-ibnd1] = atof(a[i]);
1374         AS1 = new GeomAdaptor_HSurface(GS1,UVbnd[0],UVbnd[1],UVbnd[2],UVbnd[3]);
1375         AS2 = new GeomAdaptor_HSurface(GS2,UVbnd[4],UVbnd[5],UVbnd[6],UVbnd[7]);
1376       }
1377       //
1378       if (useStart && !useBnd) {
1379         Inters.Perform(GS1,GS2,tol,UVsta[0],UVsta[1],UVsta[2],UVsta[3]);
1380       }
1381       else if (!useStart && useBnd) {
1382         Inters.Perform(AS1,AS2,tol);
1383       }
1384       else {
1385         Inters.Perform(AS1,AS2,tol,UVsta[0],UVsta[1],UVsta[2],UVsta[3]);
1386       }
1387     }//else if (n == 8 || n == 9 || n == 12 || n == 13 || n == 16 || n == 17) {
1388     else {
1389       di<<"incorrect number of arguments"<<"\n";
1390       return 1;
1391     }
1392     //
1393     if (!Inters.IsDone()) {
1394       return 1;
1395     }
1396     //
1397     char buf[1024];
1398     Standard_Integer i, aNbLines, aNbPoints; 
1399       //
1400     aNbLines = Inters.NbLines();
1401     if (aNbLines >= 2) {
1402       for (i=1; i<=aNbLines; ++i) {
1403         sprintf(buf, "%s_%d",a[1],i);
1404         Result = Inters.Line(i);
1405         const char* temp = buf; 
1406         DrawTrSurf::Set(temp,Result);
1407       }
1408     }
1409     else if (aNbLines == 1) {
1410       Result = Inters.Line(1);
1411       DrawTrSurf::Set(a[1],Result);
1412     }
1413     //
1414     aNbPoints=Inters.NbPoints();
1415     for (i=1; i<=aNbPoints; ++i) {
1416       Point=Inters.Point(i);
1417       sprintf(buf,"%s_p_%d",a[1],i);
1418       const char* temp =buf;
1419       DrawTrSurf::Set(temp, Point);
1420     }
1421   }// if (GC1.IsNull()) {
1422   //
1423   else {
1424     // Curve Surface
1425     GeomAPI_IntCS Inters(GC1,GS2);
1426     //
1427     if (!Inters.IsDone()) return 1;
1428     
1429     Standard_Integer nblines = Inters.NbSegments();
1430     Standard_Integer nbpoints = Inters.NbPoints();
1431     if ( (nblines+nbpoints) >= 2) {
1432       char newname[1024];
1433       Standard_Integer i;
1434       Standard_Integer Compt = 1;
1435       for (i = 1; i <= nblines; i++, Compt++) {
1436         sprintf(newname,"%s_%d",a[1],Compt);
1437         Result = Inters.Segment(i);
1438         const char* temp = newname; // pour portage WNT
1439         DrawTrSurf::Set(temp,Result);
1440       }
1441       for (i = 1; i <= nbpoints; i++, Compt++) {
1442         sprintf(newname,"%s_%d",a[1],i);
1443         Point = Inters.Point(i);
1444         const char* temp = newname; // pour portage WNT
1445         DrawTrSurf::Set(temp,Point);
1446       }
1447     }
1448     else if (nblines == 1) {
1449       Result = Inters.Segment(1);
1450       DrawTrSurf::Set(a[1],Result);
1451     }
1452     else if (nbpoints == 1) {
1453       Point = Inters.Point(1);
1454       DrawTrSurf::Set(a[1],Point);
1455     }
1456   }
1457
1458   dout.Flush();
1459   return 0;
1460 }
1461
1462 //=======================================================================
1463 //function : CurveCommands
1464 //purpose  : 
1465 //=======================================================================
1466 void  GeometryTest::CurveCommands(Draw_Interpretor& theCommands)
1467 {
1468   
1469   static Standard_Boolean loaded = Standard_False;
1470   if (loaded) return;
1471   loaded = Standard_True;
1472   
1473   DrawTrSurf::BasicCommands(theCommands);
1474   
1475   const char* g;
1476   
1477   g = "GEOMETRY curves creation";
1478
1479   theCommands.Add("law",
1480                   "law  name degree nbknots  knot, umult  value",
1481                   __FILE__,
1482                   polelaw,g);
1483
1484   theCommands.Add("to2d","to2d c2dname c3d [plane (XOY)]",
1485                   __FILE__,
1486                   to2d,g);
1487
1488   theCommands.Add("to3d","to3d c3dname c2d [plane (XOY)]",
1489                   __FILE__,
1490                   to3d,g);
1491
1492   theCommands.Add("gproject",
1493                   "gproject : [projectname] curve surface",
1494                   __FILE__,
1495                   gproject,g);
1496   
1497   theCommands.Add("project",
1498                   "project : no args to have help",
1499                   __FILE__,
1500                   project,g);
1501   
1502   theCommands.Add("projonplane",
1503                   "projonplane r C3d Plane [dx dy dz] [0/1]",
1504                   projonplane);
1505
1506   theCommands.Add("bisec",
1507                   "bisec result line/circle/point line/circle/point",
1508                   __FILE__,
1509                   bisec, g);
1510
1511   g = "GEOMETRY Curves and Surfaces modification";
1512
1513
1514   theCommands.Add("movelaw",
1515                   "movelaw name u  x  tx [ constraint = 0]",
1516                   __FILE__,
1517                   movelaw,g) ;
1518
1519
1520
1521   g = "GEOMETRY intersections";
1522
1523   theCommands.Add("intersect",
1524                   "intersect result surf1/curv1 surf2 [tolerance]\n\t\t  "
1525                   "intersect result surf1 surf2 [u1 v1 u2 v2] [U1F U1L V1F V1L U2F U2L V2F V2L] [tolerance]",
1526                   __FILE__,
1527                   intersection,g);
1528
1529   theCommands.Add("crvpoints",
1530                   "crvpoints result curv deflection",
1531                   __FILE__,
1532                   crvpoints,g);
1533
1534   theCommands.Add("crvtpoints",
1535                   "crvtpoints result curv deflection angular deflection - tangential deflection points",
1536                   __FILE__,
1537                   crvtpoints,g);
1538   
1539   theCommands.Add("uniformAbscissa",
1540                   "uniformAbscissa Curve nbPnt",
1541                   __FILE__,
1542                   uniformAbscissa,g);
1543
1544   theCommands.Add("uniformAbscissaEl",
1545                   "uniformAbscissaEl maxR minR nbPnt",
1546                   __FILE__,  EllipsUniformAbscissa,g);
1547
1548   theCommands.Add("mypoints",
1549                   "mypoints result curv deflection",
1550                   __FILE__,
1551                   mypoints,g);
1552   theCommands.Add("surfpoints",
1553                   "surfoints result surf deflection",
1554                   __FILE__,
1555                   surfpoints,g);
1556
1557 }
1558