src/GCPnts/GCPnts_AbscissaPoint.pxx

   1 // Created on: 1995-05-05
   2 // Created by: Modelistation
   3 // Copyright (c) 1995-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 // Dimension independant used to implement GCPnts_AbscissaPoint
  18
  19 // compute the type
  20 // and the length ratio if GCPnts_LengthParametrized
  21 #include <GCPnts_AbscissaType.hxx>
  22 #include <gp_Vec.hxx>
  23 #include <gp_Vec2d.hxx>
  24 #include <gp_Circ.hxx>
  25 #include <gp_Circ2d.hxx>
  26 #include <Precision.hxx>
  27 #include <TColStd_Array1OfReal.hxx>
  28 #include <BSplCLib.hxx>
  29
  30 static GCPnts_AbscissaType computeType( const TheCurve& C,
  31                                        Standard_Real& Ratio)
  32 {
  33   GCPnts_AbscissaType LocalType ;
  34
  35   if (C.NbIntervals(GeomAbs_CN) > 1)
  36     return GCPnts_AbsComposite;
  37
  38   switch (C.GetType()) {
  39
  40   case GeomAbs_Line:
  41     Ratio = 1.0e0 ;
  42     return GCPnts_LengthParametrized;
  43
  44   case GeomAbs_Circle:
  45     Ratio = C.Circle().Radius();
  46     return GCPnts_LengthParametrized;
  47
  48   case GeomAbs_BezierCurve:
  49     {
  50       Handle_TheBezierCurve Bz = C.Bezier();
  51       if ((Bz->NbPoles() == 2) && !(Bz->IsRational())) {
  52         Ratio = Bz->DN(0,1).Magnitude();
  53         LocalType =  GCPnts_LengthParametrized;
  54       }
  55       else
  56         LocalType = GCPnts_Parametrized;
  57       return LocalType ;
  58     }
  59   case GeomAbs_BSplineCurve:
  60     {
  61       Handle_TheBSplineCurve Bs = C.BSpline();
  62       if ((Bs->NbPoles() == 2) && !(Bs->IsRational())) {
  63         Ratio = Bs->DN(Bs->FirstParameter(),1).Magnitude();
  64         LocalType = GCPnts_LengthParametrized;
  65       }
  66       else
  67         LocalType = GCPnts_Parametrized;
  68       return LocalType ;
  69     }
  70     default:
  71       return GCPnts_Parametrized;
  72
  73   }
  74 }
  75
  76 // compute a point at distance Abscis from parameter U0
  77 // using Ui as initial guess
  78
  79 static void Compute(CPnts_AbscissaPoint& theComputer,
  80                     const TheCurve& C,
  81                     Standard_Real& Abscis,
  82                     Standard_Real& U0,
  83                     Standard_Real& Ui,
  84                     const Standard_Real EPSILON)
  85 {
  86   // test for easy solution
  87   if (Abs(Abscis) <= Precision::Confusion()) {
  88     theComputer.SetParameter(U0);
  89     return;
  90   }
  91
  92   Standard_Real Ratio = 1.;
  93   GCPnts_AbscissaType Type = computeType(C,Ratio);
  94
  95   switch (Type) {
  96   case GCPnts_LengthParametrized :
  97     theComputer.SetParameter(U0 + Abscis / Ratio);
  98     return;
  99
 100   case GCPnts_Parametrized :
 101     theComputer.Init(C);
 102     theComputer.Perform(Abscis, U0, Ui, EPSILON);
 103     return;
 104
 105   case GCPnts_AbsComposite :
 106     {
 107       Standard_Integer NbIntervals = C.NbIntervals(GeomAbs_CN);
 108       TColStd_Array1OfReal TI(1,NbIntervals+1);
 109       C.Intervals(TI,GeomAbs_CN);
 110       Standard_Real L = 0.0, sign = 1.;
 111       Standard_Integer Index = 1;
 112       BSplCLib::Hunt(TI,U0,Index);
 113       Standard_Integer Direction = 1;
 114       if (Abscis < 0) {
 115         Direction = 0;
 116         Abscis = -Abscis;
 117         sign = -1.;
 118       }
 119
 120       while ((Index >= 1) && (Index <= NbIntervals)) {
 121
 122         L = CPnts_AbscissaPoint::Length(C, U0, TI(Index+Direction));
 123         if (Abs(L - Abscis) <= Precision::Confusion()) {
 124           theComputer.SetParameter(TI(Index+Direction));
 125           return;
 126         }
 127         if(L > Abscis) {
 128           if ((Ui < TI(Index)) || (Ui > TI(Index+1))) {
 129             Ui = (Abscis / L) * (TI(Index+1) - U0);
 130             if (Direction)
 131               Ui = U0 + Ui;
 132             else
 133               Ui = U0 - Ui;
 134           }
 135           theComputer.Init(C,TI(Index),TI(Index+1));
 136           theComputer.Perform(sign*Abscis, U0, Ui, EPSILON);
 137           return;
 138         }
 139         else {
 140           U0 = TI(Index+Direction);
 141           Abscis -= L;
 142         }
 143         if (Direction)
 144           Index++;
 145         else
 146           Index--;
 147       }
 148
 149       // Push a little bit outside the limits (hairy !!!)
 150       Ui = U0 + 0.1;
 151       theComputer.Init(C,U0,U0+0.2);
 152       theComputer.Perform(sign*Abscis, U0, Ui, EPSILON);
 153       return;
 154     }
 155     break;
 156   }
 157
 158 }
 159
 160 // introduced by rbv for curvilinear parametrization
 161 // performs more apropriate tolerance managment
 162
 163 static void AdvCompute(CPnts_AbscissaPoint& theComputer,
 164                     const TheCurve& C,
 165                     Standard_Real& Abscis,
 166                     Standard_Real& U0,
 167                     Standard_Real& Ui,
 168                     const Standard_Real EPSILON)
 169 {
 170   Standard_Real Ratio = 1.0;
 171   GCPnts_AbscissaType Type = computeType(C,Ratio);
 172
 173   switch (Type) {
 174   case GCPnts_LengthParametrized :
 175     theComputer.SetParameter(U0 + Abscis / Ratio);
 176     return;
 177
 178   case GCPnts_Parametrized :
 179 //    theComputer.Init(C);
 180     theComputer.Init(C, EPSILON); //rbv's modification
 181 //
 182     theComputer.AdvPerform(Abscis, U0, Ui, EPSILON);
 183     return;
 184
 185   case GCPnts_AbsComposite :
 186     {
 187       Standard_Integer NbIntervals = C.NbIntervals(GeomAbs_CN);
 188       TColStd_Array1OfReal TI(1,NbIntervals+1);
 189       C.Intervals(TI,GeomAbs_CN);
 190       Standard_Real L = 0.0, sign = 1.;
 191       Standard_Integer Index = 1;
 192       BSplCLib::Hunt(TI,U0,Index);
 193
 194       Standard_Integer Direction = 1;
 195       if (Abscis < 0) {
 196         Direction = 0;
 197         Abscis = -Abscis;
 198         sign = -1.;
 199       }
 200
 201       if(Index == 0 && Direction > 0) {
 202         L = CPnts_AbscissaPoint::Length(C, U0, TI(Index+Direction), EPSILON);
 203         if (Abs(L - Abscis) <= /*Precision::Confusion()*/EPSILON) {
 204           theComputer.SetParameter(TI(Index+Direction));
 205           return;
 206         }
 207         if(L > Abscis) {
 208           if ( Ui > TI(Index+1) ) {
 209             Ui = (Abscis / L) * (TI(Index+1) - U0);
 210             Ui = U0 + Ui;
 211           }
 212           theComputer.Init(C,U0,TI(Index+1), EPSILON);
 213           theComputer.AdvPerform(sign*Abscis, U0, Ui, EPSILON);
 214           return;
 215         }
 216         else {
 217           U0 = TI(Index+Direction);
 218           Abscis -= L;
 219         }
 220         Index++;
 221       }
 222
 223
 224       while ((Index >= 1) && (Index <= NbIntervals)) {
 225
 226         L = CPnts_AbscissaPoint::Length(C, U0, TI(Index+Direction), EPSILON);
 227         if (Abs(L - Abscis) <= Precision::PConfusion()) {
 228           theComputer.SetParameter(TI(Index+Direction));
 229           return;
 230         }
 231         if(L > Abscis) {
 232           if ((Ui < TI(Index)) || (Ui > TI(Index+1))) {
 233             Ui = (Abscis / L) * (TI(Index+1) - U0);
 234             if (Direction)
 235               Ui = U0 + Ui;
 236             else
 237               Ui = U0 - Ui;
 238           }
 239           theComputer.Init(C,TI(Index),TI(Index+1), EPSILON);
 240           theComputer.AdvPerform(sign*Abscis, U0, Ui, EPSILON);
 241           return;
 242         }
 243         else {
 244           U0 = TI(Index+Direction);
 245           Abscis -= L;
 246         }
 247         if (Direction) {
 248           Index++;
 249
 250         }
 251         else {
 252           Index--;
 253
 254         }
 255       }
 256
 257       // Push a little bit outside the limits (hairy !!!)
 258
 259       Standard_Boolean nonperiodic = !C.IsPeriodic();
 260       Ui = U0 + sign*0.1;
 261       Standard_Real U1 = U0 + sign*.2;
 262       if(nonperiodic) {
 263         if(sign > 0) {
 264           Ui = Min(Ui,C.LastParameter());
 265           U1 = Min(U1, C.LastParameter());
 266         }
 267         else {
 268           Ui = Max(Ui,C.FirstParameter());
 269           U1 = Max(U1, C.FirstParameter());
 270         }
 271       }
 272
 273       theComputer.Init(C, U0, U1, EPSILON);
 274       theComputer.AdvPerform(sign*Abscis, U0, Ui, EPSILON);
 275       return;
 276     }
 277     break;
 278   }
 279
 280 }
 281
 282 //=======================================================================
 283 //function : Length
 284 //purpose  :
 285 //=======================================================================
 286
 287 Standard_Real GCPnts_AbscissaPoint::Length(const TheCurve& C)
 288 {
 289   return GCPnts_AbscissaPoint::Length(C,C.FirstParameter(),
 290                                       C.LastParameter());
 291 }
 292
 293 //=======================================================================
 294 //function : Length
 295 //purpose  :
 296 //=======================================================================
 297
 298 Standard_Real GCPnts_AbscissaPoint::Length(const TheCurve& C,
 299                                            const Standard_Real Tol)
 300 {
 301   return GCPnts_AbscissaPoint::Length(C,C.FirstParameter(),
 302                                       C.LastParameter(),Tol);
 303 }
 304
 305
 306 //=======================================================================
 307 //function : Length
 308 //purpose  :
 309 //=======================================================================
 310
 311 Standard_Real GCPnts_AbscissaPoint::Length(const TheCurve& C,
 312                                            const Standard_Real U1,
 313                                            const Standard_Real U2)
 314 {
 315   Standard_Real Ratio = 1.0;
 316   GCPnts_AbscissaType Type = computeType(C,Ratio);
 317   switch (Type) {
 318
 319   case GCPnts_LengthParametrized:
 320     return Abs(U2-U1) * Ratio;
 321
 322   case GCPnts_Parametrized:
 323     return CPnts_AbscissaPoint::Length(C, U1, U2);
 324
 325   case GCPnts_AbsComposite:
 326     {
 327       Standard_Integer NbIntervals = C.NbIntervals(GeomAbs_CN);
 328       TColStd_Array1OfReal TI(1,NbIntervals+1);
 329       C.Intervals(TI,GeomAbs_CN);
 330       Standard_Real UU1 = Min(U1, U2);
 331       Standard_Real UU2 = Max(U1, U2);
 332       Standard_Real L = 0.0;
 333       for(Standard_Integer Index = 1; Index <= NbIntervals; Index++) {
 334         if (TI(Index)   > UU2) break;
 335         if (TI(Index+1) < UU1) continue;
 336         L += CPnts_AbscissaPoint::Length(C,
 337                                          Max(TI(Index),UU1),
 338                                          Min(TI(Index+1),UU2));
 339       }
 340       return L;
 341     }
 342   }
 343   return RealLast();
 344 }
 345
 346 //=======================================================================
 347 //function : Length
 348 //purpose  :
 349 //=======================================================================
 350
 351 Standard_Real GCPnts_AbscissaPoint::Length(const TheCurve& C,
 352                                            const Standard_Real U1,
 353                                            const Standard_Real U2,
 354                                            const Standard_Real Tol)
 355 {
 356   Standard_Real Ratio = 1.0;
 357   GCPnts_AbscissaType Type = computeType(C,Ratio);
 358   switch (Type) {
 359
 360   case GCPnts_LengthParametrized:
 361     return Abs(U2-U1) * Ratio;
 362
 363   case GCPnts_Parametrized:
 364     return CPnts_AbscissaPoint::Length(C, U1, U2, Tol);
 365
 366   case GCPnts_AbsComposite:
 367     {
 368       Standard_Integer NbIntervals = C.NbIntervals(GeomAbs_CN);
 369       TColStd_Array1OfReal TI(1,NbIntervals+1);
 370       C.Intervals(TI,GeomAbs_CN);
 371       Standard_Real UU1 = Min(U1, U2);
 372       Standard_Real UU2 = Max(U1, U2);
 373       Standard_Real L = 0.0;
 374       for(Standard_Integer Index = 1; Index <= NbIntervals; Index++) {
 375         if (TI(Index)   > UU2) break;
 376         if (TI(Index+1) < UU1) continue;
 377         L += CPnts_AbscissaPoint::Length(C,
 378                                          Max(TI(Index),UU1),
 379                                          Min(TI(Index+1),UU2),
 380                                          Tol);
 381       }
 382       return L;
 383     }
 384   }
 385   return RealLast();
 386 }
 387
 388
 389 //=======================================================================
 390 //function : GCPnts_AbscissaPoint
 391 //purpose  :
 392 //=======================================================================
 393
 394 GCPnts_AbscissaPoint::GCPnts_AbscissaPoint
 395                                   (const TheCurve& C,
 396                                    const Standard_Real   Abscissa,
 397                                    const Standard_Real   U0)
 398 {
 399   Standard_Real L = GCPnts_AbscissaPoint::Length(C);
 400   if (L < Precision::Confusion()) {
 401     throw Standard_ConstructionError();
 402   }
 403   Standard_Real Abscis = Abscissa;
 404   Standard_Real UU0 = U0;
 405   Standard_Real UUi = U0 +
 406     (Abscis / L) * (C.LastParameter() - C.FirstParameter());
 407   Compute(myComputer, C, Abscis, UU0, UUi,
 408           C.Resolution(Precision::Confusion()));
 409 }
 410
 411 //=======================================================================
 412 //function : GCPnts_AbscissaPoint
 413 //purpose  : rbv for curvilinear parametrization
 414 //=======================================================================
 415
 416 GCPnts_AbscissaPoint::GCPnts_AbscissaPoint
 417                                   (const Standard_Real Tol,
 418                                    const TheCurve& C,
 419                                    const Standard_Real   Abscissa,
 420                                    const Standard_Real   U0)
 421 {
 422   Standard_Real L = GCPnts_AbscissaPoint::Length(C, Tol);
 423 /*  if (L < Precision::Confusion()) {
 424     cout<<"FirstParameter = "<<C.FirstParameter()<<endl;
 425     cout<<"LastParameter = "<<C.LastParameter()<<endl;
 426     throw Standard_ConstructionError("GCPnts_AbscissaPoint::GCPnts_AbscissaPoint");
 427   }
 428 */
 429   Standard_Real Abscis = Abscissa;
 430   Standard_Real UU0 = U0;
 431   Standard_Real UUi;
 432   if (L >= Precision::Confusion())
 433     UUi= U0 +
 434       (Abscis / L) * (C.LastParameter() - C.FirstParameter());
 435   else UUi = U0;
 436
 437   AdvCompute(myComputer, C, Abscis, UU0, UUi, Tol);
 438 }
 439
 440 //=======================================================================
 441 //function : GCPnts_AbscissaPoint
 442 //purpose  :
 443 //=======================================================================
 444
 445 GCPnts_AbscissaPoint::GCPnts_AbscissaPoint
 446                                   (const TheCurve& C,
 447                                    const Standard_Real   Abscissa,
 448                                    const Standard_Real   U0,
 449                                    const Standard_Real   Ui)
 450 {
 451   Standard_Real Abscis = Abscissa;
 452   Standard_Real UU0 = U0;
 453   Standard_Real UUi = Ui;
 454   Compute(myComputer, C, Abscis, UU0, UUi,
 455           C.Resolution(Precision::Confusion()));
 456 }
 457
 458 //=======================================================================
 459 //function : GCPnts_AbscissaPoint
 460 //purpose  : rbv for curvilinear parametrization
 461 //=======================================================================
 462
 463 GCPnts_AbscissaPoint::GCPnts_AbscissaPoint
 464                                   (const TheCurve& C,
 465                                    const Standard_Real   Abscissa,
 466                                    const Standard_Real   U0,
 467                                    const Standard_Real   Ui,
 468                                    const Standard_Real   Tol)
 469 {
 470   Standard_Real Abscis = Abscissa;
 471   Standard_Real UU0 = U0;
 472   Standard_Real UUi = Ui;
 473   AdvCompute(myComputer, C, Abscis, UU0, UUi, Tol);
 474 }