1// Copyright (c) 1995-1999 Matra Datavision 2// Copyright (c) 1999-2014 OPEN CASCADE SAS 3// 4// This file is part of Open CASCADE Technology software library. 5// 6// This library is free software; you can redistribute it and/or modify it under 7// the terms of the GNU Lesser General Public License version 2.1 as published 8// by the Free Software Foundation, with special exception defined in the file 9// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT 10// distribution for complete text of the license and disclaimer of any warranty. 11// 12// Alternatively, this file may be used under the terms of Open CASCADE 13// commercial license or contractual agreement. 14 15#include <Approx_ParametrizationType.hxx> 16#include Approx_ParLeastSquareOfMyGradient_hxx 17#include <TColStd_Array1OfReal.hxx> 18#include <TColgp_Array1OfPnt.hxx> 19#include <TColgp_Array1OfPnt2d.hxx> 20#include <gp_Pnt.hxx> 21#include <gp_Pnt2d.hxx> 22#include <gp_Vec.hxx> 23#include <gp_Vec2d.hxx> 24#include <TColgp_Array1OfVec.hxx> 25#include <TColgp_Array1OfVec2d.hxx> 26#include <AppParCurves_Constraint.hxx> 27#include <AppParCurves_HArray1OfConstraintCouple.hxx> 28#include <AppParCurves_MultiPoint.hxx> 29#include <Precision.hxx> 30#include <math_IntegerVector.hxx> 31#include <math_Gauss.hxx> 32#include <math_Uzawa.hxx> 33#include <Approx_MCurvesToBSpCurve.hxx> 34#include <AppParCurves_ConstraintCouple.hxx> 35 36#include <stdio.h> 37 38#ifdef OCCT_DEBUG 39static Standard_Boolean mydebug = Standard_False; 40 41#include <Geom_BezierCurve.hxx> 42#include <Geom2d_BezierCurve.hxx> 43#ifdef DRAW 44#include <DrawTrSurf.hxx> 45#include <Draw.hxx> 46#include <Draw_Appli.hxx> 47#endif 48 49static void DUMP(const MultiLine& Line) 50{ 51 Standard_Integer i, j, nbP2d, nbP3d, firstP, lastP; 52 gp_Pnt P1; 53 gp_Pnt2d P12d; 54 55 firstP = LineTool::FirstPoint(Line); 56 lastP = LineTool::LastPoint(Line); 57 58 nbP3d = LineTool::NbP3d(Line); 59 nbP2d = LineTool::NbP2d(Line); 60 Standard_Integer mynbP3d=nbP3d, mynbP2d=nbP2d; 61 if (nbP3d == 0) mynbP3d = 1; 62 if (nbP2d == 0) mynbP2d = 1; 63 64 TColgp_Array1OfPnt tabP(1, mynbP3d); 65 TColgp_Array1OfPnt2d tabP2d(1, mynbP2d); 66 67 std::cout <<"DUMP de la MultiLine entre "<<firstP <<" et "<<lastP<<": "<<std::endl; 68 for (i = firstP; i <= lastP; i++) { 69 if (nbP3d != 0 && nbP2d != 0) LineTool::Value(Line, i, tabP, tabP2d); 70 else if (nbP2d != 0) LineTool::Value(Line, i, tabP2d); 71 else if (nbP3d != 0) LineTool::Value(Line, i, tabP); 72 73 std::cout << "point "<<i<<":"<< std::endl; 74 for (j = 1; j <= nbP3d; j++) { 75 P1 = tabP(j); 76 std::cout <<P1.X()<<" "<<P1.Y()<<" "<<P1.Z()<<std::endl; 77 } 78 for (j = 1; j <= nbP2d; j++) { 79 P12d = tabP2d(j); 80 std::cout <<P12d.X()<<" "<<P12d.Y()<<std::endl; 81 } 82 } 83 84} 85 86static void DUMP(const AppParCurves_MultiCurve& C) { 87 static Standard_Integer nbappel = 0; 88 Standard_Integer i; 89 Standard_Integer nbpoles = C.NbPoles(); 90 91 Handle(Geom_BezierCurve) BSp; 92 Handle(Geom2d_BezierCurve) BSp2d; 93 94 TColgp_Array1OfPnt tabPoles(1, nbpoles); 95 TColgp_Array1OfPnt2d tabPoles2d(1, nbpoles); 96 char solname[100]; 97 98 nbappel++; 99 for (i = 1; i <= C.NbCurves(); i++) { 100 if (C.Dimension(i) == 3) { 101 C.Curve(i, tabPoles); 102 BSp = new Geom_BezierCurve(tabPoles); 103 sprintf(solname, "%s%i%s_%i", "c", i, "3d", nbappel); 104#ifdef DRAW 105 char* Temp = solname; 106 DrawTrSurf::Set(Temp, BSp); 107// DrawTrSurf::Set(solname, BSp); 108#endif 109 } 110 else { 111 C.Curve(i, tabPoles2d); 112 BSp2d = new Geom2d_BezierCurve(tabPoles2d); 113 sprintf(solname, "%s%i%s_%i", "c", i, "2d", nbappel); 114#ifdef DRAW 115 char* Temp = solname; 116 DrawTrSurf::Set(Temp, BSp2d); 117// DrawTrSurf::Set(solname, BSp2d); 118#endif 119 } 120 } 121#ifdef DRAW 122 dout.Flush(); 123#endif 124} 125 126#endif 127 128static Standard_Boolean CheckMultiCurve(const AppParCurves_MultiCurve& theMultiCurve, 129 const MultiLine& theLine, 130 const Standard_Integer theIndfirst, 131 const Standard_Integer theIndlast, 132 Standard_Integer& theIndbad) 133{ 134 const Standard_Integer nbp3d = LineTool::NbP3d(theLine); 135 const Standard_Integer nbp2d = LineTool::NbP2d(theLine); 136 137 const Standard_Real coeff = 4.; //2*2 138 139 if (nbp3d > 1) //only simple cases are analysed 140 return Standard_True; 141 142 const Standard_Real MinScalProd = -0.9; 143 const Standard_Real SqTol3d = Precision::SquareConfusion(); 144 145 theIndbad = 0; 146 Standard_Integer indbads [4]; 147 indbads[1] = indbads[2] = indbads[3] = 0; 148 149 Standard_Integer NbCur = theMultiCurve.NbCurves(); 150 Standard_Boolean LoopFound = Standard_False; 151 152 Standard_Integer aNbP3d = Max(nbp3d, 1); 153 Standard_Integer aNbP2d = Max(nbp2d, 1); 154 155 TColgp_Array1OfPnt tabP(1, aNbP3d); 156 TColgp_Array1OfPnt2d tabP2d(1, aNbP2d); 157 158#ifdef DRAW 159 char* name = new char[100]; 160 Standard_Integer nbbc = 1; 161 Standard_Integer indc = 1; 162#endif 163 if (theMultiCurve.Dimension(1) == 3 /*myNbP3d == 1*/) 164 { 165 TColgp_Array1OfPnt aPoles(1, theMultiCurve.NbPoles()); 166 theMultiCurve.Curve(1, aPoles); 167#ifdef DRAW 168 Handle(Geom_Curve) theBezier = new Geom_BezierCurve(aPoles); 169 sprintf(name, "bc3d_%d_%d", indc, nbbc); 170 DrawTrSurf::Set(name, theBezier); 171#endif 172 gp_Vec FirstVec, SecondVec; 173 Standard_Integer indp = 2; 174 while (indp <= aPoles.Upper()) 175 { 176 FirstVec = gp_Vec(aPoles(1), aPoles(indp++)); 177 Standard_Real aLength = FirstVec.Magnitude(); 178 if (aLength > gp::Resolution()) 179 { 180 FirstVec /= aLength; 181 break; 182 } 183 } 184 gp_Pnt MidPnt = aPoles(indp-1); 185 //for (Standard_Integer k = 3; k <= aPoles.Upper(); k++) 186 while (indp <= aPoles.Upper()) 187 { 188 SecondVec = gp_Vec(MidPnt, aPoles(indp)); 189 Standard_Real aLength = SecondVec.Magnitude(); 190 if (aLength <= gp::Resolution()) 191 { 192 indp++; 193 continue; 194 } 195 SecondVec /= aLength; 196 Standard_Real ScalProd = FirstVec * SecondVec; 197 if (ScalProd < MinScalProd) 198 { 199#ifdef DRAW 200 std::cout<<"ScalProd("<<indp-2<<","<<indp-1<<")-("<<indp-1<<","<<indp<<") = "<<ScalProd<<std::endl; 201#endif 202 LoopFound = Standard_True; 203 break; 204 } 205 FirstVec = SecondVec; 206 MidPnt = aPoles(indp); 207 indp++; 208 } 209 //Check: may be it is a real loop 210 if (LoopFound) 211 { 212#ifdef DRAW 213 for (Standard_Integer ipoint = theIndfirst; ipoint <= theIndlast; ipoint++) 214 { 215 LineTool::Value(theLine, ipoint, tabP); 216 gp_Pnt aPnt = tabP(1); 217 sprintf(name, "p%d", ipoint); 218 DrawTrSurf::Set(name, aPnt); 219 } 220#endif 221 for (Standard_Integer FirstInd = theIndfirst; 222 FirstInd <= theIndlast - 2; FirstInd++) 223 { 224 LineTool::Value(theLine, FirstInd, tabP); 225 gp_Pnt FirstPnt = tabP(1); 226 for (Standard_Integer k = FirstInd+1; k < theIndlast; k++) 227 { 228 LineTool::Value(theLine, k, tabP); 229 gp_Pnt Pnt1 = tabP(1); 230 LineTool::Value(theLine, k+1, tabP); 231 gp_Pnt Pnt2 = tabP(1); 232 if (FirstPnt.SquareDistance(Pnt1) <= SqTol3d || 233 FirstPnt.SquareDistance(Pnt2) <= SqTol3d) 234 { 235 LoopFound = Standard_False; 236 break; 237 } 238 gp_Vec Vec1(FirstPnt, Pnt1); 239 Vec1.Normalize(); 240 gp_Vec Vec2(FirstPnt, Pnt2); 241 Vec2.Normalize(); 242 Standard_Real ScalProd = Vec1 * Vec2; 243 if (ScalProd < MinScalProd) 244 { 245 LoopFound = Standard_False; 246 break; 247 } 248 } 249 if (LoopFound == Standard_False) 250 break; 251 } 252 } 253 if (LoopFound) 254 { 255 //search <indbad> 256 Standard_Real MaxSqDist = 0.; 257 Standard_Real MinSqDist = RealLast(); 258 for (Standard_Integer k = theIndfirst+1; k <= theIndlast; k++) 259 { 260 LineTool::Value(theLine, k-1, tabP); 261 gp_Pnt PrevPnt = tabP(1); 262 LineTool::Value(theLine, k, tabP); 263 gp_Pnt CurPnt = tabP(1); 264 Standard_Real aSqDist = PrevPnt.SquareDistance(CurPnt); 265 if (aSqDist > MaxSqDist) 266 { 267 MaxSqDist = aSqDist; 268 indbads[1] = k; 269 } 270 if (aSqDist > gp::Resolution() && 271 aSqDist < MinSqDist) 272 MinSqDist = aSqDist; 273 } 274 Standard_Real Relation = MaxSqDist / MinSqDist; 275 if (Relation < coeff) 276 LoopFound = Standard_False; 277 else 278 for (Standard_Integer indcur = 2; indcur <= NbCur; indcur++) 279 { 280 MaxSqDist = 0.; 281 for (Standard_Integer k = theIndfirst+1; k <= theIndlast; k++) 282 { 283 LineTool::Value(theLine, k-1, tabP2d); 284 gp_Pnt2d PrevPnt = tabP2d(indcur-1); 285 LineTool::Value(theLine, k, tabP2d); 286 gp_Pnt2d CurPnt = tabP2d(indcur-1); 287 Standard_Real aSqDist = PrevPnt.SquareDistance(CurPnt); 288 if (aSqDist > MaxSqDist) 289 { 290 MaxSqDist = aSqDist; 291 indbads[indcur] = k; 292 } 293 } 294 } 295 } 296 } //if (myNbP3d == 1) 297 else //2d case 298 { 299 TColgp_Array1OfPnt2d aPoles2d(1, theMultiCurve.NbPoles()); 300 theMultiCurve.Curve(1, aPoles2d); 301#ifdef DRAW 302 Handle(Geom2d_Curve) theBezier2d = new Geom2d_BezierCurve(aPoles2d); 303 sprintf(name, "bc2d_%d_%d", indc, nbbc); 304 DrawTrSurf::Set(name, theBezier2d); 305#endif 306 const Standard_Real aSqNormToler = Epsilon(1.0)*Epsilon(1.0); 307 gp_Vec2d FirstVec(aPoles2d(1), aPoles2d(2)), SecondVec; 308 Standard_Real aVecSqNorm = FirstVec.SquareMagnitude(); 309 if (aVecSqNorm < aSqNormToler) 310 { 311 theIndbad = theIndfirst + 1; 312 return Standard_False; 313 } 314 315 FirstVec /= Sqrt(aSqNormToler); 316 gp_Pnt2d MidPnt = aPoles2d(2); 317 for (Standard_Integer k = 3; k <= aPoles2d.Upper(); k++) 318 { 319 SecondVec.SetXY(aPoles2d(k).XY() - MidPnt.XY()); 320 aVecSqNorm = SecondVec.SquareMagnitude(); 321 if (aVecSqNorm < aSqNormToler) 322 { 323 theIndbad = theIndfirst + k - 1; 324 return Standard_False; 325 } 326 327 SecondVec /= Sqrt(aVecSqNorm); 328 Standard_Real ScalProd = FirstVec * SecondVec; 329 if (ScalProd < MinScalProd) 330 { 331#ifdef DRAW 332 std::cout<<"ScalProd("<<k-2<<","<<k-1<<")-("<<k-1<<","<<k<<") = "<<ScalProd<<std::endl; 333#endif 334 LoopFound = Standard_True; 335 break; 336 } 337 FirstVec = SecondVec; 338 MidPnt = aPoles2d(k); 339 } 340 //Check: may be it is a real loop 341 if (LoopFound) 342 { 343#ifdef DRAW 344 for (Standard_Integer ipoint = theIndfirst; ipoint <= theIndlast; ipoint++) 345 { 346 LineTool::Value(theLine, ipoint, tabP2d); 347 gp_Pnt2d aPnt2d = tabP2d(1); 348 sprintf(name, "p%d", ipoint); 349 DrawTrSurf::Set(name, aPnt2d); 350 } 351#endif 352 for (Standard_Integer FirstInd = theIndfirst; 353 FirstInd <= theIndlast - 2; FirstInd++) 354 { 355 LineTool::Value(theLine, FirstInd, tabP2d); 356 gp_Pnt2d FirstPnt = tabP2d(1); 357 for (Standard_Integer k = FirstInd+1; k < theIndlast; k++) 358 { 359 LineTool::Value(theLine, k, tabP2d); 360 gp_Pnt2d Pnt1 = tabP2d(1); 361 LineTool::Value(theLine, k+1, tabP2d); 362 gp_Pnt2d Pnt2 = tabP2d(1); 363 if (FirstPnt.SquareDistance(Pnt1) <= SqTol3d || 364 FirstPnt.SquareDistance(Pnt2) <= SqTol3d) 365 { 366 LoopFound = Standard_False; 367 break; 368 } 369 gp_Vec2d Vec1(FirstPnt, Pnt1); 370 Vec1.Normalize(); 371 gp_Vec2d Vec2(FirstPnt, Pnt2); 372 Vec2.Normalize(); 373 Standard_Real ScalProd = Vec1 * Vec2; 374 if (ScalProd < MinScalProd) 375 { 376 LoopFound = Standard_False; 377 break; 378 } 379 } 380 if (LoopFound == Standard_False) 381 break; 382 } 383 } 384 if (LoopFound) 385 { 386 //search <indbad> 387 for (Standard_Integer indcur = 1; indcur <= NbCur; indcur++) 388 { 389 Standard_Real MaxSqDist = 0.; 390 Standard_Real MinSqDist = RealLast(); 391 for (Standard_Integer k = theIndfirst+1; k <= theIndlast; k++) 392 { 393 LineTool::Value(theLine, k-1, tabP2d); 394 gp_Pnt2d PrevPnt = tabP2d(indcur); 395 LineTool::Value(theLine, k, tabP2d); 396 gp_Pnt2d CurPnt = tabP2d(indcur); 397 Standard_Real aSqDist = PrevPnt.SquareDistance(CurPnt); 398 if (aSqDist > MaxSqDist) 399 { 400 MaxSqDist = aSqDist; 401 indbads[indcur] = k; 402 } 403 if (aSqDist > gp::Resolution() && 404 aSqDist < MinSqDist) 405 MinSqDist = aSqDist; 406 } 407 Standard_Real Relation = MaxSqDist / MinSqDist; 408 if (Relation < coeff) 409 LoopFound = Standard_False; 410 } 411 } 412 } 413 414 //Define <indbad> 415 for (Standard_Integer i = 1; i <= 3; i++) 416 if (indbads[i] != 0) 417 { 418 theIndbad = indbads[i]; 419 break; 420 } 421 422 if (!LoopFound) 423 theIndbad = 0; 424 425 return (!LoopFound); 426} 427 428void Approx_ComputeLine::FirstTangencyVector(const MultiLine& Line, 429 const Standard_Integer index, 430 math_Vector& V) const 431{ 432 433 Standard_Integer i, j, nbP2d, nbP3d; 434 nbP3d = LineTool::NbP3d(Line); 435 nbP2d = LineTool::NbP2d(Line); 436 Standard_Integer mynbP3d=nbP3d, mynbP2d=nbP2d; 437 if (nbP3d == 0) mynbP3d = 1; 438 if (nbP2d == 0) mynbP2d = 1; 439 Standard_Boolean Ok=Standard_False; 440 TColgp_Array1OfVec TabV(1, mynbP3d); 441 TColgp_Array1OfVec2d TabV2d(1, mynbP2d); 442 443 if (nbP3d != 0 && nbP2d != 0) 444 Ok = LineTool::Tangency(Line, index, TabV, TabV2d); 445 else if (nbP2d != 0) 446 Ok = LineTool::Tangency(Line, index, TabV2d); 447 else if (nbP3d != 0) 448 Ok = LineTool::Tangency(Line, index, TabV); 449 450 if (Ok) { 451 if (nbP3d != 0) { 452 j = 1; 453 for (i = TabV.Lower(); i <= TabV.Upper(); i++) { 454 V(j) = TabV(i).X(); 455 V(j+1) = TabV(i).Y(); 456 V(j+2) = TabV(i).Z(); 457 j += 3; 458 } 459 } 460 if (nbP2d != 0) { 461 j = nbP3d*3+1; 462 for (i = TabV2d.Lower(); i <= TabV2d.Upper(); i++) { 463 V(j) = TabV2d(i).X(); 464 V(j+1) = TabV2d(i).Y(); 465 j += 2; 466 } 467 } 468 } 469 else { 470 471 // recherche d un vecteur tangent par construction d une parabole: 472 AppParCurves_Constraint firstC, lastC; 473 firstC = lastC = AppParCurves_PassPoint; 474 Standard_Integer nbpoles = 3; 475 math_Vector mypar(index, index+2); 476 Parameters(Line, index, index+2, mypar); 477 Approx_ParLeastSquareOfMyGradient 478 LSQ(Line, index, index+2, firstC, lastC, mypar, nbpoles); 479 AppParCurves_MultiCurve C = LSQ.BezierValue(); 480 481 gp_Pnt myP; 482 gp_Vec myV; 483 gp_Pnt2d myP2d; 484 gp_Vec2d myV2d; 485 j = 1; 486 for (i = 1; i <= nbP3d; i++) { 487 C.D1(i, 0.0, myP, myV); 488 V(j) = myV.X(); 489 V(j+1) = myV.Y(); 490 V(j+2) = myV.Z(); 491 j += 3; 492 } 493 j = nbP3d*3+1; 494 for (i = nbP3d+1; i <= nbP3d+nbP2d; i++) { 495 C.D1(i, 0.0, myP2d, myV2d); 496 V(j) = myV2d.X(); 497 V(j+1) = myV2d.Y(); 498 j += 2; 499 } 500 } 501} 502 503 504void Approx_ComputeLine::LastTangencyVector(const MultiLine& Line, 505 const Standard_Integer index, 506 math_Vector& V) const 507{ 508 Standard_Integer i, j, nbP2d, nbP3d; 509 nbP3d = LineTool::NbP3d(Line); 510 nbP2d = LineTool::NbP2d(Line); 511 Standard_Integer mynbP3d=nbP3d, mynbP2d=nbP2d; 512 if (nbP3d == 0) mynbP3d = 1; 513 if (nbP2d == 0) mynbP2d = 1; 514 Standard_Boolean Ok=Standard_False; 515 TColgp_Array1OfVec TabV(1, mynbP3d); 516 TColgp_Array1OfVec2d TabV2d(1, mynbP2d); 517 518 519 if (nbP3d != 0 && nbP2d != 0) 520 Ok = LineTool::Tangency(Line, index, TabV, TabV2d); 521 else if (nbP2d != 0) 522 Ok = LineTool::Tangency(Line, index, TabV2d); 523 else if (nbP3d != 0) 524 Ok = LineTool::Tangency(Line, index, TabV); 525 526 if (Ok) { 527 if (nbP3d != 0) { 528 j = 1; 529 for (i = TabV.Lower(); i <= TabV.Upper(); i++) { 530 V(j) = TabV(i).X(); 531 V(j+1) = TabV(i).Y(); 532 V(j+2) = TabV(i).Z(); 533 j += 3; 534 } 535 } 536 if (nbP2d != 0) { 537 j = nbP3d*3+1; 538 for (i = TabV2d.Lower(); i <= TabV2d.Upper(); i++) { 539 V(j) = TabV2d(i).X(); 540 V(j+1) = TabV2d(i).Y(); 541 j += 2; 542 } 543 } 544 } 545 else { 546 547 // recherche d un vecteur tangent par construction d une parabole: 548 AppParCurves_Constraint firstC, lastC; 549 firstC = lastC = AppParCurves_PassPoint; 550 Standard_Integer nbpoles = 3; 551 math_Vector mypar(index-2, index); 552 Parameters(Line, index-2, index, mypar); 553 Approx_ParLeastSquareOfMyGradient 554 LSQ(Line, index-2, index, firstC, lastC, mypar, nbpoles); 555 AppParCurves_MultiCurve C = LSQ.BezierValue(); 556 557 gp_Pnt myP; 558 gp_Vec myV; 559 gp_Pnt2d myP2d; 560 gp_Vec2d myV2d; 561 j = 1; 562 for (i = 1; i <= nbP3d; i++) { 563 C.D1(i, 1.0, myP, myV); 564 V(j) = myV.X(); 565 V(j+1) = myV.Y(); 566 V(j+2) = myV.Z(); 567 j += 3; 568 } 569 j = nbP3d*3+1; 570 for (i = nbP3d+1; i <= nbP3d+nbP2d; i++) { 571 C.D1(i, 1.0, myP2d, myV2d); 572 V(j) = myV2d.X(); 573 V(j+1) = myV2d.Y(); 574 j += 2; 575 } 576 } 577 578} 579 580 581 582Standard_Real Approx_ComputeLine:: 583 SearchFirstLambda(const MultiLine& Line, 584 const math_Vector& TheParam, 585 const math_Vector& V, 586 const Standard_Integer index) const 587{ 588 589 // dq/dw = lambda* V = (p2-p1)/(u2-u1) 590 591 Standard_Integer nbP2d, nbP3d; 592 gp_Pnt P1, P2; 593 gp_Pnt2d P12d, P22d; 594 nbP3d = LineTool::NbP3d(Line); 595 nbP2d = LineTool::NbP2d(Line); 596 Standard_Integer mynbP3d=nbP3d, mynbP2d=nbP2d; 597 if (nbP3d == 0) mynbP3d = 1; 598 if (nbP2d == 0) mynbP2d = 1; 599 TColgp_Array1OfPnt tabP1(1, mynbP3d), tabP2(1, mynbP3d); 600 TColgp_Array1OfPnt2d tabP12d(1, mynbP2d), tabP22d(1, mynbP2d); 601 602 603 if (nbP3d != 0 && nbP2d != 0) LineTool::Value(Line, index, tabP1, tabP12d); 604 else if (nbP2d != 0) LineTool::Value(Line, index, tabP12d); 605 else if (nbP3d != 0) LineTool::Value(Line, index, tabP1); 606 607 if (nbP3d != 0 && nbP2d != 0) LineTool::Value(Line, index+1, tabP2, tabP22d); 608 else if (nbP2d != 0) LineTool::Value(Line, index+1, tabP22d); 609 else if (nbP3d != 0) LineTool::Value(Line, index+1, tabP2); 610 611 612 Standard_Real U1 = TheParam(index), U2 = TheParam(index+1); 613 Standard_Real lambda, S; 614 Standard_Integer low = V.Lower(); 615 616 if (nbP3d != 0) { 617 P1 = tabP1(1); 618 P2 = tabP2(1); 619 gp_Vec P1P2(P1, P2), myV; 620 myV.SetCoord(V(low), V(low+1), V(low+2)); 621 lambda = (P1P2.Magnitude())/(myV.Magnitude()*(U2-U1)); 622 S = (P1P2.Dot(myV)> 0.0) ? 1.0 : -1.0; 623 } 624 else { 625 P12d = tabP12d(1); 626 P22d = tabP22d(1); 627 gp_Vec2d P1P2(P12d, P22d), myV; 628 myV.SetCoord(V(low), V(low+1)); 629 lambda = (P1P2.Magnitude())/(myV.Magnitude()*(U2-U1)); 630 S = (P1P2.Dot(myV)> 0.0) ? 1.0 : -1.0; 631 } 632 return (S*lambda); 633 634} 635 636 637Standard_Real Approx_ComputeLine:: 638 SearchLastLambda(const MultiLine& Line, 639 const math_Vector& TheParam, 640 const math_Vector& V, 641 const Standard_Integer index) const 642{ 643 // dq/dw = lambda* V = (p2-p1)/(u2-u1) 644 645 Standard_Integer nbP2d, nbP3d; 646 gp_Pnt P1, P2; 647 gp_Pnt2d P12d, P22d; 648 nbP3d = LineTool::NbP3d(Line); 649 nbP2d = LineTool::NbP2d(Line); 650 Standard_Integer mynbP3d=nbP3d, mynbP2d=nbP2d; 651 if (nbP3d == 0) mynbP3d = 1; 652 if (nbP2d == 0) mynbP2d = 1; 653 TColgp_Array1OfPnt tabP(1, mynbP3d), tabP2(1, mynbP3d); 654 TColgp_Array1OfPnt2d tabP2d(1, mynbP2d), tabP22d(1, mynbP2d); 655 656 if (nbP3d != 0 && nbP2d != 0) LineTool::Value(Line, index-1, tabP, tabP2d); 657 else if (nbP2d != 0) LineTool::Value(Line, index-1, tabP2d); 658 else if (nbP3d != 0) LineTool::Value(Line, index-1, tabP); 659 660 if (nbP3d != 0 && nbP2d != 0) LineTool::Value(Line, index, tabP2, tabP22d); 661 else if (nbP2d != 0) LineTool::Value(Line, index, tabP22d); 662 else if (nbP3d != 0) LineTool::Value(Line, index, tabP2); 663 664 665 Standard_Real U1 = TheParam(index-1), U2 = TheParam(index); 666 Standard_Real lambda, S; 667 Standard_Integer low = V.Lower(); 668 669 if (nbP3d != 0) { 670 P1 = tabP(1); 671 P2 = tabP2(1); 672 gp_Vec P1P2(P1, P2), myV; 673 myV.SetCoord(V(low), V(low+1), V(low+2)); 674 lambda = (P1P2.Magnitude())/(myV.Magnitude()*(U2-U1)); 675 S = (P1P2.Dot(myV)> 0.0) ? 1.0 : -1.0; 676 } 677 else { 678 P12d = tabP2d(1); 679 P22d = tabP22d(1); 680 gp_Vec2d P1P2(P12d, P22d), myV; 681 myV.SetCoord(V(low), V(low+1)); 682 lambda = (P1P2.Magnitude())/(myV.Magnitude()*(U2-U1)); 683 S = (P1P2.Dot(myV)> 0.0) ? 1.0 : -1.0; 684 } 685 686 return (S*lambda); 687} 688 689 690 691Approx_ComputeLine::Approx_ComputeLine 692 (const MultiLine& Line, 693 const math_Vector& Parameters, 694 const Standard_Integer degreemin, 695 const Standard_Integer degreemax, 696 const Standard_Real Tolerance3d, 697 const Standard_Real Tolerance2d, 698 const Standard_Integer NbIterations, 699 const Standard_Boolean cutting, 700 const Standard_Boolean Squares) 701: myMultiLineNb (0), 702 myIsClear (Standard_False) 703{ 704 myfirstParam = new TColStd_HArray1OfReal(Parameters.Lower(), 705 Parameters.Upper()); 706 for (Standard_Integer i = Parameters.Lower(); i <= Parameters.Upper(); i++) { 707 myfirstParam->SetValue(i, Parameters(i)); 708 } 709 myConstraints = new AppParCurves_HArray1OfConstraintCouple(1, 2); 710 Par = Approx_IsoParametric; 711 mydegremin = degreemin; 712 mydegremax = degreemax; 713 mytol3d = Tolerance3d; 714 mytol2d = Tolerance2d; 715 mysquares = Squares; 716 mycut = cutting; 717 myitermax = NbIterations; 718 alldone = Standard_False; 719 myfirstC = AppParCurves_TangencyPoint; 720 mylastC = AppParCurves_TangencyPoint; 721 Perform(Line); 722} 723 724 725Approx_ComputeLine::Approx_ComputeLine 726 (const math_Vector& Parameters, 727 const Standard_Integer degreemin, 728 const Standard_Integer degreemax, 729 const Standard_Real Tolerance3d, 730 const Standard_Real Tolerance2d, 731 const Standard_Integer NbIterations, 732 const Standard_Boolean cutting, 733 const Standard_Boolean Squares) 734: myMultiLineNb (0), 735 myIsClear (Standard_False) 736{ 737 myfirstParam = new TColStd_HArray1OfReal(Parameters.Lower(), 738 Parameters.Upper()); 739 for (Standard_Integer i = Parameters.Lower(); i <= Parameters.Upper(); i++) { 740 myfirstParam->SetValue(i, Parameters(i)); 741 } 742 myfirstC = AppParCurves_TangencyPoint; 743 mylastC = AppParCurves_TangencyPoint; 744 myConstraints = new AppParCurves_HArray1OfConstraintCouple(1, 2); 745 Par = Approx_IsoParametric; 746 mydegremin = degreemin; 747 mydegremax = degreemax; 748 mytol3d = Tolerance3d; 749 mytol2d = Tolerance2d; 750 mysquares = Squares; 751 mycut = cutting; 752 myitermax = NbIterations; 753 alldone = Standard_False; 754} 755 756Approx_ComputeLine::Approx_ComputeLine 757 (const Standard_Integer degreemin, 758 const Standard_Integer degreemax, 759 const Standard_Real Tolerance3d, 760 const Standard_Real Tolerance2d, 761 const Standard_Integer NbIterations, 762 const Standard_Boolean cutting, 763 const Approx_ParametrizationType parametrization, 764 const Standard_Boolean Squares) 765: myMultiLineNb (0), 766 myIsClear (Standard_False) 767{ 768 myConstraints = new AppParCurves_HArray1OfConstraintCouple(1, 2); 769 Par = parametrization; 770 mydegremin = degreemin; 771 mydegremax = degreemax; 772 mytol3d = Tolerance3d; 773 mytol2d = Tolerance2d; 774 mysquares = Squares; 775 mycut = cutting; 776 myitermax = NbIterations; 777 myfirstC = AppParCurves_TangencyPoint; 778 mylastC = AppParCurves_TangencyPoint; 779 alldone = Standard_False; 780} 781 782 783Approx_ComputeLine::Approx_ComputeLine 784 (const MultiLine& Line, 785 const Standard_Integer degreemin, 786 const Standard_Integer degreemax, 787 const Standard_Real Tolerance3d, 788 const Standard_Real Tolerance2d, 789 const Standard_Integer NbIterations, 790 const Standard_Boolean cutting, 791 const Approx_ParametrizationType parametrization, 792 const Standard_Boolean Squares) 793: myMultiLineNb (0), 794 myIsClear (Standard_False) 795{ 796 myConstraints = new AppParCurves_HArray1OfConstraintCouple(1, 2); 797 alldone = Standard_False; 798 mydegremin = degreemin; 799 mydegremax = degreemax; 800 mytol3d = Tolerance3d; 801 mytol2d = Tolerance2d; 802 mysquares = Squares; 803 mycut = cutting; 804 myitermax = NbIterations; 805 Par = parametrization; 806 myfirstC = AppParCurves_TangencyPoint; 807 mylastC = AppParCurves_TangencyPoint; 808 809 Perform(Line); 810} 811 812 813 814void Approx_ComputeLine::Perform(const MultiLine& Line) 815{ 816#ifdef OCCT_DEBUG 817 if (mydebug) DUMP(Line); 818#endif 819 if (!myIsClear) 820 { 821 myMultiCurves.Clear(); 822 myPar.Clear(); 823 Tolers3d.Clear(); 824 Tolers2d.Clear(); 825 myMultiLineNb = 0; 826 } 827 else myIsClear = Standard_False; 828 829 Standard_Integer i, nbp, Thefirstpt, Thelastpt, oldlastpt; 830 Standard_Boolean Finish = Standard_False, 831 begin = Standard_True, Ok = Standard_False, 832 GoUp = Standard_False, Interpol; 833 Standard_Real thetol3d, thetol2d; 834 Approx_Status MyStatus; 835// gp_Vec V13d, V23d; 836// gp_Vec2d V2d; 837 Thefirstpt = LineTool::FirstPoint(Line); 838 Thelastpt = LineTool::LastPoint(Line); 839 Standard_Integer myfirstpt = Thefirstpt; 840 Standard_Integer mylastpt = Thelastpt; 841 842 AppParCurves_ConstraintCouple myCouple1(myfirstpt, myfirstC); 843 AppParCurves_ConstraintCouple myCouple2(mylastpt, mylastC); 844 myConstraints->SetValue(1, myCouple1); 845 myConstraints->SetValue(2, myCouple2); 846 847 math_Vector TheParam(Thefirstpt, Thelastpt); 848 849 850 if (!mycut) { 851 if(myfirstParam.IsNull()) { 852 Parameters(Line, Thefirstpt, Thelastpt, TheParam); 853 } 854 else { 855 for (i = myfirstParam->Lower(); i <= myfirstParam->Upper(); i++) { 856 TheParam(i+Thefirstpt-1) = myfirstParam->Value(i); 857 } 858 } 859 TheMultiCurve = AppParCurves_MultiCurve(); 860 MultiLine anOtherLine0; 861 Standard_Boolean isOtherLine0Made = Standard_False; 862 Standard_Integer indbad = 0; 863 alldone = Compute(Line, myfirstpt, mylastpt, TheParam, thetol3d, thetol2d, indbad); 864 if (indbad != 0) 865 { 866 isOtherLine0Made = LineTool::MakeMLOneMorePoint(Line, myfirstpt, mylastpt, indbad, anOtherLine0); 867 } 868 if (isOtherLine0Made) 869 { 870 myIsClear = Standard_True; 871 //++myMultiLineNb; 872 Perform(anOtherLine0); 873 alldone = Standard_True; 874 } 875 if(!alldone && TheMultiCurve.NbCurves() > 0) { 876#ifdef OCCT_DEBUG 877 if (mydebug) DUMP(TheMultiCurve); 878#endif 879 myMultiCurves.Append(TheMultiCurve); 880 Tolers3d.Append(currenttol3d); 881 Tolers2d.Append(currenttol2d); 882 Standard_Integer mylen = mylastpt-myfirstpt+1; 883 Standard_Integer myParLen = myParameters->Length(); 884 Standard_Integer aLen = (myParLen > mylen)? myParLen : mylen; 885 Handle(TColStd_HArray1OfReal) ThePar = 886 new TColStd_HArray1OfReal(myfirstpt, myfirstpt+aLen-1); 887 for (i = 0; i < aLen; i++) 888 ThePar->SetValue(myfirstpt+i, myParameters->Value(myParameters->Lower()+i)); 889 myPar.Append(ThePar); 890 } 891 } 892 else { 893 while (!Finish) { 894 oldlastpt = mylastpt; 895 // Gestion du decoupage de la multiline pour approximer: 896 if(!begin) { 897 if (!GoUp) { 898 if (Ok) { 899 // Calcul de la partie a approximer. 900 myfirstpt = mylastpt; 901 mylastpt = Thelastpt; 902 if (myfirstpt == Thelastpt) { 903 Finish = Standard_True; 904 alldone = Standard_True; 905 return; 906 } 907 } 908 else { 909 nbp = mylastpt - myfirstpt + 1; 910 MyStatus = LineTool::WhatStatus(Line, myfirstpt, mylastpt); 911 if (MyStatus == Approx_NoPointsAdded && nbp <= mydegremax+1) { 912 Interpol = ComputeCurve(Line, myfirstpt, mylastpt); 913 if (Interpol) { 914 if (mylastpt == Thelastpt) { 915 Finish = Standard_True; 916 alldone = Standard_True; 917 return; 918 } 919 } 920 } 921 mylastpt = Standard_Integer((myfirstpt + mylastpt)/2); 922 } 923 } 924 GoUp = Standard_False; 925 } 926 927 // Verification du nombre de points restants par rapport au degre 928 // demande. 929 // ============================================================== 930 nbp = mylastpt - myfirstpt + 1; 931 MyStatus = LineTool::WhatStatus(Line, myfirstpt, mylastpt); 932 if (nbp <= mydegremax+5 ) { 933 // Rajout necessaire de points si possible. 934 // ======================================== 935 GoUp = Standard_False; 936 Ok = Standard_True; 937 if (MyStatus == Approx_PointsAdded) { 938 // Appel recursif du decoupage: 939 GoUp = Standard_True; 940 941 MultiLine anOtherLine1 = LineTool::MakeMLBetween(Line, myfirstpt, mylastpt, nbp-1); 942 943 Standard_Integer nbpdsotherligne = LineTool::FirstPoint (anOtherLine1) - LineTool::LastPoint (anOtherLine1); 944 945 //-- Si MakeML a echoue on retourne une ligne vide 946 if ((nbpdsotherligne == 0) || myMultiLineNb >= 3) 947 { 948 //-- cout<<" ** ApproxComputeLine MakeML Echec ** LBR lbr "<<endl; 949 if (myfirstpt == mylastpt) break; // Pour etre sur de ne pas 950 // planter la station !! 951 myCouple1.SetIndex(myfirstpt); 952 myCouple2.SetIndex(mylastpt); 953 myConstraints->SetValue(1, myCouple1); 954 myConstraints->SetValue(2, myCouple2); 955 956 math_Vector Param(myfirstpt, mylastpt); 957 Approx_ParametrizationType SavePar = Par; 958 Par = Approx_IsoParametric; 959 Parameters(Line, myfirstpt, mylastpt, Param); 960 TheMultiCurve = AppParCurves_MultiCurve(); 961 MultiLine anOtherLine2; 962 Standard_Boolean isOtherLine2Made = Standard_False; 963 Standard_Integer indbad = 0; 964 Ok = Compute(Line, myfirstpt, mylastpt, Param, thetol3d, thetol2d, indbad); 965 if (indbad != 0) 966 { 967 isOtherLine2Made = LineTool::MakeMLOneMorePoint(Line, myfirstpt, mylastpt, indbad, anOtherLine2); 968 } 969 if (isOtherLine2Made) 970 { 971 myIsClear = Standard_True; 972 //++myMultiLineNb; 973 Par = SavePar; 974 Perform(anOtherLine2); 975 Ok = Standard_True; 976 } 977 978 if (!Ok) { 979 Standard_Real tt3d = currenttol3d, tt2d = currenttol2d; 980 Handle(TColStd_HArray1OfReal) saveParameters = myParameters; 981 AppParCurves_MultiCurve saveMultiCurve = TheMultiCurve; 982 983 if(SavePar != Approx_IsoParametric) 984 Par = SavePar; 985 else 986 Par = Approx_ChordLength; 987 988 Parameters(Line, myfirstpt, mylastpt, Param); 989 isOtherLine2Made = Standard_False; 990 indbad = 0; 991 Ok = Compute(Line, myfirstpt, mylastpt, Param, thetol3d, thetol2d, indbad); 992 if (indbad != 0) 993 { 994 isOtherLine2Made = LineTool::MakeMLOneMorePoint (Line, myfirstpt, mylastpt, indbad, anOtherLine2); 995 } 996 if (isOtherLine2Made) 997 { 998 myIsClear = Standard_True; 999 //++myMultiLineNb; 1000 Perform (anOtherLine2); 1001 Ok = Standard_True; 1002 } 1003 1004 if (!Ok && tt3d <= currenttol3d && tt2d <= currenttol2d) { 1005 currenttol3d = tt3d; currenttol2d = tt2d; 1006 myParameters = saveParameters; 1007 TheMultiCurve = saveMultiCurve; 1008 } 1009 } 1010 Par = SavePar; 1011 if (myfirstpt == Thelastpt) 1012 { 1013 Finish = Standard_True; 1014 alldone = Standard_True; 1015 return; 1016 } 1017 1018 oldlastpt = mylastpt; 1019 if (!Ok) { 1020 tolreached = Standard_False; 1021 if (TheMultiCurve.NbCurves() == 0) { 1022 myMultiCurves.Clear(); 1023 return; 1024 } 1025#ifdef OCCT_DEBUG 1026 if (mydebug) DUMP(TheMultiCurve); 1027#endif 1028 MultiLine anOtherLine3; 1029 Standard_Boolean isOtherLine3Made = Standard_False; 1030 Standard_Integer indbad2 = 0; 1031 if (!CheckMultiCurve(TheMultiCurve, Line, 1032 myfirstpt, mylastpt, 1033 indbad2)) 1034 { 1035 isOtherLine3Made = LineTool::MakeMLOneMorePoint (Line, myfirstpt, mylastpt, indbad2, anOtherLine3); 1036 } 1037 if (isOtherLine3Made) 1038 { 1039 myIsClear = Standard_True; 1040 //++myMultiLineNb; 1041 Perform(anOtherLine3); 1042 myfirstpt = mylastpt; 1043 mylastpt = Thelastpt; 1044 } 1045 else 1046 { 1047 myMultiCurves.Append(TheMultiCurve); 1048 Tolers3d.Append(currenttol3d); 1049 Tolers2d.Append(currenttol2d); 1050 Standard_Integer mylen = oldlastpt-myfirstpt+1; 1051 Standard_Integer myParLen = myParameters->Length(); 1052 Standard_Integer aLen = (myParLen > mylen)? myParLen : mylen; 1053 Handle(TColStd_HArray1OfReal) ThePar = 1054 new TColStd_HArray1OfReal(myfirstpt, myfirstpt+aLen-1); 1055 for (i = 0; i < aLen; i++) 1056 ThePar->SetValue(myfirstpt+i, myParameters->Value(myParameters->Lower()+i)); 1057 myPar.Append(ThePar); 1058 } 1059 } 1060 myfirstpt = oldlastpt; 1061 mylastpt = Thelastpt; 1062 1063 } 1064 else 1065 { 1066 myIsClear = Standard_True; 1067 ++myMultiLineNb; 1068 Perform(anOtherLine1); 1069 myfirstpt = mylastpt; 1070 mylastpt = Thelastpt; 1071 } 1072 } 1073 1074 if (MyStatus == Approx_NoPointsAdded && !begin) { 1075 // On rend la meilleure approximation obtenue precedemment. 1076 // ======================================================== 1077 GoUp = Standard_True; 1078 tolreached = Standard_False; 1079 if (TheMultiCurve.NbCurves() == 0) { 1080 myMultiCurves.Clear(); 1081 return; 1082 } 1083#ifdef OCCT_DEBUG 1084 if (mydebug) DUMP(TheMultiCurve); 1085#endif 1086 myMultiCurves.Append(TheMultiCurve); 1087 Tolers3d.Append(currenttol3d); 1088 Tolers2d.Append(currenttol2d); 1089 Standard_Integer mylen = oldlastpt-myfirstpt+1; 1090 Standard_Integer myParLen = myParameters->Length(); 1091 Standard_Integer aLen = (myParLen > mylen)? myParLen : mylen; 1092 Handle(TColStd_HArray1OfReal) ThePar = 1093 new TColStd_HArray1OfReal(myfirstpt, myfirstpt+aLen-1); 1094 for (i = 0; i < aLen; i++) 1095 ThePar->SetValue(myfirstpt+i, myParameters->Value(myParameters->Lower()+i)); 1096 myPar.Append(ThePar); 1097 1098 myfirstpt = oldlastpt; 1099 mylastpt = Thelastpt; 1100 } 1101 1102 else if (MyStatus == Approx_NoApproximation) { 1103 // On ne fait pas d approximation entre myfirstpt et mylastpt. 1104 // =========================================================== 1105 // On stocke pour pouvoir en informer l utilisateur. 1106 GoUp = Standard_True; 1107 myfirstpt = mylastpt; 1108 mylastpt = Thelastpt; 1109 } 1110 } 1111 1112 if (myfirstpt == Thelastpt) { 1113 Finish = Standard_True; 1114 alldone = Standard_True; 1115 return; 1116 } 1117 if (!GoUp) { 1118 if (myfirstpt == mylastpt) break; // Pour etre sur de ne pas 1119 // planter la station !! 1120 myCouple1.SetIndex(myfirstpt); 1121 myCouple2.SetIndex(mylastpt); 1122 myConstraints->SetValue(1, myCouple1); 1123 myConstraints->SetValue(2, myCouple2); 1124 1125 // Calcul des parametres sur ce nouvel intervalle. 1126 // On recupere les parametres initiaux lors du decoupage. 1127 1128 math_Vector Param(myfirstpt, mylastpt); 1129 if (begin) { 1130 if(myfirstParam.IsNull()) { 1131 Parameters(Line, myfirstpt, mylastpt, Param); 1132 } 1133 else { 1134 for (i = myfirstParam->Lower(); i <= myfirstParam->Upper(); i++) { 1135 Param(i) = myfirstParam->Value(i); 1136 } 1137 myfirstParam.Nullify(); 1138 } 1139 TheParam = Param; 1140 begin = Standard_False; 1141 } 1142 else { 1143 Standard_Real pfirst = TheParam.Value(myfirstpt); 1144 Standard_Real plast = TheParam.Value(mylastpt); 1145 for (i = myfirstpt; i <= mylastpt; i++) { 1146 Param(i) = (TheParam.Value(i)-pfirst)/(plast-pfirst); 1147 } 1148 } 1149 1150 TheMultiCurve = AppParCurves_MultiCurve(); 1151 Standard_Integer indbad = 0; 1152 Ok = Compute(Line, myfirstpt, mylastpt, Param, thetol3d, thetol2d, indbad); 1153 if (myfirstpt == Thelastpt) 1154 { 1155 Finish = Standard_True; 1156 alldone = Standard_True; 1157 return; 1158 } 1159 } 1160 } 1161 } 1162} 1163 1164 1165 1166const TColStd_Array1OfReal& Approx_ComputeLine::Parameters(const Standard_Integer Index) const 1167{ 1168 return (myPar.Value(Index))->Array1(); 1169} 1170 1171 1172Standard_Integer Approx_ComputeLine::NbMultiCurves()const 1173{ 1174 return myMultiCurves.Length(); 1175} 1176 1177AppParCurves_MultiCurve& Approx_ComputeLine::ChangeValue(const Standard_Integer Index) 1178{ 1179 return myMultiCurves.ChangeValue(Index); 1180} 1181 1182 1183const AppParCurves_MultiCurve& Approx_ComputeLine::Value(const Standard_Integer Index) 1184const 1185{ 1186 return myMultiCurves.Value(Index); 1187} 1188 1189 1190const AppParCurves_MultiBSpCurve& Approx_ComputeLine::SplineValue() 1191{ 1192 Approx_MCurvesToBSpCurve Trans; 1193 Trans.Perform(myMultiCurves); 1194 myspline = Trans.Value(); 1195 return myspline; 1196} 1197 1198 1199 1200 1201 1202void Approx_ComputeLine::Parameters(const MultiLine& Line, 1203 const Standard_Integer firstP, 1204 const Standard_Integer lastP, 1205 math_Vector& TheParameters) const 1206{ 1207 Standard_Integer i, j, nbP2d, nbP3d; 1208 Standard_Real dist; 1209 1210 if(Par == Approx_ChordLength || Par == Approx_Centripetal) 1211 { 1212 nbP3d = LineTool::NbP3d(Line); 1213 nbP2d = LineTool::NbP2d(Line); 1214 Standard_Integer mynbP3d = nbP3d, mynbP2d = nbP2d; 1215 if(nbP3d == 0) mynbP3d = 1; 1216 if(nbP2d == 0) mynbP2d = 1; 1217 1218 TheParameters(firstP) = 0.0; 1219 dist = 0.0; 1220 TColgp_Array1OfPnt tabP(1, mynbP3d); 1221 TColgp_Array1OfPnt tabPP(1, mynbP3d); 1222 TColgp_Array1OfPnt2d tabP2d(1, mynbP2d); 1223 TColgp_Array1OfPnt2d tabPP2d(1, mynbP2d); 1224 1225 for(i = firstP + 1; i <= lastP; i++) 1226 { 1227 if(nbP3d != 0 && nbP2d != 0) LineTool::Value(Line, i - 1, tabP, tabP2d); 1228 else if(nbP2d != 0) LineTool::Value(Line, i - 1, tabP2d); 1229 else if(nbP3d != 0) LineTool::Value(Line, i - 1, tabP); 1230 1231 if(nbP3d != 0 && nbP2d != 0) LineTool::Value(Line, i, tabPP, tabPP2d); 1232 else if(nbP2d != 0) LineTool::Value(Line, i, tabPP2d); 1233 else if(nbP3d != 0) LineTool::Value(Line, i, tabPP); 1234 dist = 0; 1235 for(j = 1; j <= nbP3d; j++) 1236 { 1237 const gp_Pnt &aP1 = tabP(j), 1238 &aP2 = tabPP(j); 1239 dist += aP2.SquareDistance(aP1); 1240 } 1241 for(j = 1; j <= nbP2d; j++) 1242 { 1243 const gp_Pnt2d &aP12d = tabP2d(j), 1244 &aP22d = tabPP2d(j); 1245 1246 dist += aP22d.SquareDistance(aP12d); 1247 } 1248 1249 dist = Sqrt(dist); 1250 if(Par == Approx_ChordLength) 1251 { 1252 TheParameters(i) = TheParameters(i - 1) + dist; 1253 } 1254 else 1255 {// Par == Approx_Centripetal 1256 TheParameters(i) = TheParameters(i - 1) + Sqrt(dist); 1257 } 1258 } 1259 for(i = firstP; i <= lastP; i++) TheParameters(i) /= TheParameters(lastP); 1260 } 1261 else { 1262 for (i = firstP; i <= lastP; i++) { 1263 TheParameters(i) = (Standard_Real(i)-firstP)/ 1264 (Standard_Real(lastP)-Standard_Real(firstP)); 1265 } 1266 } 1267} 1268 1269 1270Standard_Boolean Approx_ComputeLine::Compute(const MultiLine& Line, 1271 const Standard_Integer fpt, 1272 const Standard_Integer lpt, 1273 math_Vector& Para, 1274 Standard_Real& TheTol3d, 1275 Standard_Real& TheTol2d, 1276 Standard_Integer& indbad) 1277{ 1278 indbad = 0; 1279 Standard_Integer deg, i; 1280 Standard_Boolean mydone; 1281 Standard_Real Fv; 1282 Standard_Integer nbp = lpt-fpt+1; 1283 1284 math_Vector ParSav(Para.Lower(), Para.Upper()); 1285 for (i = Para.Lower(); i <= Para.Upper(); i++) { 1286 ParSav(i) = Para(i); 1287 } 1288 Standard_Integer Mdegmax = mydegremax; 1289 if(nbp < Mdegmax+5 && mycut) { 1290 Mdegmax = nbp - 5; 1291 } 1292 if(Mdegmax < mydegremin) { 1293 Mdegmax = mydegremin; 1294 } 1295 1296 currenttol3d = currenttol2d = RealLast(); 1297 for (deg = Min(nbp-1,mydegremin); deg <= Mdegmax; deg++) { 1298 AppParCurves_MultiCurve mySCU(deg+1); 1299 if (mysquares) { 1300 Approx_ParLeastSquareOfMyGradient SQ(Line, fpt, lpt, 1301 myfirstC, mylastC, Para, deg+1); 1302 mydone = SQ.IsDone(); 1303 mySCU = SQ.BezierValue(); 1304 SQ.Error(Fv, TheTol3d, TheTol2d); 1305 } 1306 else { 1307 Approx_MyGradient GRAD(Line, fpt, lpt, myConstraints, 1308 Para, deg, mytol3d, mytol2d, myitermax); 1309 mydone = GRAD.IsDone(); 1310 mySCU = GRAD.Value(); 1311 if (mySCU.NbCurves() == 0) 1312 continue; 1313 TheTol3d = GRAD.MaxError3d(); 1314 TheTol2d = GRAD.MaxError2d(); 1315 } 1316 Standard_Real uu1 = Para(Para.Lower()), uu2; 1317 Standard_Boolean restau = Standard_False; 1318 for ( i = Para.Lower()+1; i <= Para.Upper(); i++) { 1319 uu2 = Para(i); 1320 if (uu2 <= uu1) { 1321 restau = Standard_True; 1322// cout << "restau = Standard_True" << endl; 1323 break; 1324 } 1325 uu1 = uu2; 1326 } 1327 if (restau) { 1328 for (i = Para.Lower(); i <= Para.Upper(); i++) { 1329 Para(i) = ParSav(i); 1330 } 1331 } 1332 if (mydone) { 1333 if (TheTol3d <= mytol3d && TheTol2d <= mytol2d) { 1334 // Stockage de la multicurve approximee. 1335 tolreached = Standard_True; 1336#ifdef OCCT_DEBUG 1337 if (mydebug) DUMP(mySCU); 1338#endif 1339 if (!CheckMultiCurve(mySCU, Line, 1340 fpt, lpt, 1341 indbad)) 1342 { 1343 return Standard_False; 1344 } 1345 else 1346 { 1347 myMultiCurves.Append(mySCU); 1348 // Stockage des parametres de la partie de MultiLine approximee: 1349 // A ameliorer !! (bq trop de recopies) 1350 Handle(TColStd_HArray1OfReal) ThePar = 1351 new TColStd_HArray1OfReal(Para.Lower(), Para.Upper()); 1352 for (i = Para.Lower(); i <= Para.Upper(); i++) { 1353 ThePar->SetValue(i, Para(i)); 1354 } 1355 myPar.Append(ThePar); 1356 Tolers3d.Append(TheTol3d); 1357 Tolers2d.Append(TheTol2d); 1358 return Standard_True; 1359 } 1360 } 1361 } 1362 1363 if (TheTol3d <= currenttol3d && TheTol2d <= currenttol2d) { 1364 TheMultiCurve = mySCU; 1365 currenttol3d = TheTol3d; 1366 currenttol2d = TheTol2d; 1367 myParameters = new TColStd_HArray1OfReal(Para.Lower(), Para.Upper()); 1368 for (i = Para.Lower(); i <= Para.Upper(); i++) { 1369 myParameters->SetValue(i, Para(i)); 1370 } 1371 } 1372 1373 } 1374 1375 return Standard_False; 1376} 1377 1378 1379 1380 1381Standard_Boolean Approx_ComputeLine::ComputeCurve(const MultiLine& Line, 1382 const Standard_Integer firstpt, 1383 const Standard_Integer lastpt) 1384{ 1385 Standard_Integer i, j, nbP3d, nbP2d, deg; 1386 gp_Vec V13d, V23d; 1387 gp_Vec2d V12d, V22d; 1388 gp_Pnt P1, P2; 1389 gp_Pnt2d P12d, P22d; 1390 Standard_Boolean Tangent1, Tangent2, mydone= Standard_False; 1391#ifdef OCCT_DEBUG 1392 Standard_Boolean Parallel; 1393#endif 1394 Standard_Integer myfirstpt = firstpt, mylastpt = lastpt; 1395 Standard_Integer nbp = lastpt-firstpt+1, Kopt = 0; 1396 math_Vector Para(firstpt, lastpt); 1397 1398 Parameters(Line, firstpt, lastpt, Para); 1399 1400 nbP3d = LineTool::NbP3d(Line); 1401 nbP2d = LineTool::NbP2d(Line); 1402 Standard_Integer mynbP3d = nbP3d, mynbP2d = nbP2d; 1403 if (nbP3d == 0) mynbP3d = 1 ; 1404 if (nbP2d == 0) mynbP2d = 1 ; 1405 1406 1407 TColgp_Array1OfVec tabV1(1, mynbP3d), tabV2(1, mynbP3d); 1408 TColgp_Array1OfPnt tabP1(1, mynbP3d), tabP2(1, mynbP3d), tabP(1, mynbP3d); 1409 TColgp_Array1OfVec2d tabV12d(1, mynbP2d), tabV22d(1, mynbP2d); 1410 TColgp_Array1OfPnt2d tabP12d(1, mynbP2d), tabP22d(1, mynbP2d), tabP2d(1, mynbP2d); 1411 1412 if (nbP3d != 0 && nbP2d != 0) { 1413 LineTool::Value(Line, myfirstpt,tabP1,tabP12d); 1414 LineTool::Value(Line, mylastpt,tabP2,tabP22d); 1415 Tangent1 = LineTool::Tangency(Line, myfirstpt,tabV1,tabV12d); 1416 Tangent2 = LineTool::Tangency(Line, mylastpt,tabV2,tabV22d); 1417 } 1418 else if (nbP2d != 0) { 1419 LineTool::Value(Line, myfirstpt,tabP12d); 1420 LineTool::Value(Line, mylastpt,tabP22d); 1421 Tangent1 = LineTool::Tangency(Line, myfirstpt, tabV12d); 1422 Tangent2 = LineTool::Tangency(Line, mylastpt, tabV22d); 1423 } 1424 else { 1425 LineTool::Value(Line, myfirstpt,tabP1); 1426 LineTool::Value(Line, mylastpt,tabP2); 1427 Tangent1 = LineTool::Tangency(Line, myfirstpt, tabV1); 1428 Tangent2 = LineTool::Tangency(Line, mylastpt, tabV2); 1429 } 1430 1431 if (Tangent1) Kopt++; 1432 if (Tangent2) Kopt++; 1433 1434 1435 if (nbp == 2) { 1436 // S il n y a que 2 points, on verifie quand meme que les tangentes sont 1437 // alignees. 1438#ifdef OCCT_DEBUG 1439 Parallel = Standard_True; 1440#endif 1441 if (Tangent1) { 1442 for (i = 1; i <= nbP3d; i++) { 1443 gp_Vec PVec(tabP1(i), tabP2(i)); 1444 V13d = tabV1(i); 1445 if (!PVec.IsParallel(V13d, Precision::Angular())) { 1446#ifdef OCCT_DEBUG 1447 Parallel = Standard_False; 1448#endif 1449 break; 1450 } 1451 } 1452 for (i = 1; i <= nbP2d; i++) { 1453 gp_Vec2d PVec2d(tabP12d(i), tabP22d(i)); 1454 V12d = tabV12d(i); 1455 if (!PVec2d.IsParallel(V12d, Precision::Angular())) { 1456#ifdef OCCT_DEBUG 1457 Parallel = Standard_False; 1458#endif 1459 break; 1460 } 1461 } 1462 } 1463 1464 if (Tangent2) { 1465 for (i = 1; i <= nbP3d; i++) { 1466 gp_Vec PVec(tabP1(i), tabP2(i)); 1467 V23d = tabV2(i); 1468 if (!PVec.IsParallel(V23d, Precision::Angular())) { 1469#ifdef OCCT_DEBUG 1470 Parallel = Standard_False; 1471#endif 1472 break; 1473 } 1474 } 1475 for (i = 1; i <= nbP2d; i++) { 1476 gp_Vec2d PVec2d(tabP12d(i), tabP22d(i)); 1477 V22d = tabV22d(i); 1478 if (!PVec2d.IsParallel(V22d, Precision::Angular())) { 1479#ifdef OCCT_DEBUG 1480 Parallel = Standard_False; 1481#endif 1482 break; 1483 } 1484 } 1485 } 1486 1487#ifdef OCCT_DEBUG 1488 if (!Parallel) { 1489 if (mydebug) std::cout <<"droite mais tangentes pas vraiment paralleles!!"<< std::endl; 1490 } 1491#endif 1492 AppParCurves_MultiCurve mySCU(mydegremin+1); 1493 if (nbP3d != 0 && nbP2d != 0) { 1494 AppParCurves_MultiPoint MPole1(tabP1, tabP12d); 1495 AppParCurves_MultiPoint MPole2(tabP2, tabP22d); 1496 mySCU.SetValue(1, MPole1); 1497 mySCU.SetValue(mydegremin+1, MPole2); 1498 for (i = 2; i <= mydegremin; i++) { 1499 for (j = 1; j<= nbP3d; j++) { 1500 P1 = tabP1(j); 1501 P2 = tabP2(j); 1502 tabP(j).SetXYZ(P1.XYZ()+(i-1)*(P2.XYZ()-P1.XYZ())/mydegremin); 1503 } 1504 for (j = 1; j<= nbP2d; j++) { 1505 P12d = tabP12d(j); 1506 P22d = tabP22d(j); 1507 tabP2d(j).SetXY(P12d.XY()+(i-1)*(P22d.XY()-P12d.XY())/mydegremin); 1508 } 1509 AppParCurves_MultiPoint MPole(tabP, tabP2d); 1510 mySCU.SetValue(i, MPole); 1511 } 1512 1513 } 1514 else if (nbP3d != 0) { 1515 AppParCurves_MultiPoint MPole1(tabP1); 1516 AppParCurves_MultiPoint MPole2(tabP2); 1517 mySCU.SetValue(1, MPole1); 1518 mySCU.SetValue(mydegremin+1, MPole2); 1519 for (i = 2; i <= mydegremin; i++) { 1520 for (j = 1; j<= nbP3d; j++) { 1521 P1 = tabP1(j); 1522 P2 = tabP2(j); 1523 tabP(j).SetXYZ(P1.XYZ()+(i-1)*(P2.XYZ()-P1.XYZ())/mydegremin); 1524 } 1525 AppParCurves_MultiPoint MPole(tabP); 1526 mySCU.SetValue(i, MPole); 1527 } 1528 } 1529 else if (nbP2d != 0) { 1530 AppParCurves_MultiPoint MPole1(tabP12d); 1531 AppParCurves_MultiPoint MPole2(tabP22d); 1532 mySCU.SetValue(1, MPole1); 1533 mySCU.SetValue(mydegremin+1, MPole2); 1534 for (i = 2; i <= mydegremin; i++) { 1535 for (j = 1; j<= nbP2d; j++) { 1536 P12d = tabP12d(j); 1537 P22d = tabP22d(j); 1538 tabP2d(j).SetXY(P12d.XY()+(i-1)*(P22d.XY()-P12d.XY())/mydegremin); 1539 } 1540 AppParCurves_MultiPoint MPole(tabP2d); 1541 mySCU.SetValue(i, MPole); 1542 } 1543 } 1544 mydone = Standard_True; 1545 // Stockage de la multicurve approximee. 1546 tolreached = Standard_True; 1547#ifdef OCCT_DEBUG 1548 if (mydebug) DUMP(mySCU); 1549#endif 1550 myMultiCurves.Append(mySCU); 1551 Handle(TColStd_HArray1OfReal) ThePar = new TColStd_HArray1OfReal(Para.Lower(), Para.Upper()); 1552 for (i = Para.Lower(); i <= Para.Upper(); i++) { 1553 ThePar->SetValue(i, Para(i)); 1554 } 1555 myPar.Append(ThePar); 1556 Tolers3d.Append(Precision::Confusion()); 1557 Tolers2d.Append(Precision::PConfusion()); 1558 return mydone; 1559 } 1560 1561 // avec les tangentes. 1562 deg = nbp+1; 1563 AppParCurves_MultiCurve mySCU(deg+1); 1564 AppParCurves_Constraint Cons = AppParCurves_TangencyPoint; 1565 Standard_Real lambda1, lambda2; 1566 math_Vector V1(1, nbP3d*3+nbP2d*2); 1567 math_Vector V2(1, nbP3d*3+nbP2d*2); 1568 FirstTangencyVector(Line, myfirstpt, V1); 1569 lambda1 = SearchFirstLambda(Line, Para, V1, myfirstpt); 1570 1571 LastTangencyVector(Line, mylastpt, V2); 1572 lambda2 = SearchLastLambda(Line, Para, V2, mylastpt); 1573 1574 Approx_ParLeastSquareOfMyGradient 1575 LSQ(Line, myfirstpt, mylastpt, 1576 Cons, Cons, Para, deg+1); 1577 1578 lambda1 = lambda1/deg; 1579 lambda2 = lambda2/deg; 1580 LSQ.Perform(Para, V1, V2, lambda1, lambda2); 1581 mydone = LSQ.IsDone(); 1582 mySCU = LSQ.BezierValue(); 1583 1584 if (mydone) { 1585 Standard_Real Fv, TheTol3d, TheTol2d; 1586 LSQ.Error(Fv, TheTol3d, TheTol2d); 1587 1588 // Stockage de la multicurve approximee. 1589 tolreached = Standard_True; 1590#ifdef OCCT_DEBUG 1591 if (mydebug) DUMP(mySCU); 1592#endif 1593 myMultiCurves.Append(mySCU); 1594 Handle(TColStd_HArray1OfReal) ThePar = 1595 new TColStd_HArray1OfReal(Para.Lower(), Para.Upper()); 1596 for (i = Para.Lower(); i <= Para.Upper(); i++) { 1597 ThePar->SetValue(i, Para(i)); 1598 } 1599 myPar.Append(ThePar); 1600 Tolers3d.Append(TheTol3d); 1601 Tolers2d.Append(TheTol2d); 1602 return Standard_True; 1603 } 1604 return mydone; 1605} 1606 1607 1608 1609void Approx_ComputeLine::Init(const Standard_Integer degreemin, 1610 const Standard_Integer degreemax, 1611 const Standard_Real Tolerance3d, 1612 const Standard_Real Tolerance2d, 1613 const Standard_Integer NbIterations, 1614 const Standard_Boolean cutting, 1615 const Approx_ParametrizationType parametrization, 1616 const Standard_Boolean Squares) 1617{ 1618 mydegremin = degreemin; 1619 mydegremax = degreemax; 1620 mytol3d = Tolerance3d; 1621 mytol2d = Tolerance2d; 1622 Par = parametrization; 1623 mysquares = Squares; 1624 mycut = cutting; 1625 myitermax = NbIterations; 1626} 1627 1628 1629 1630void Approx_ComputeLine::SetDegrees(const Standard_Integer degreemin, 1631 const Standard_Integer degreemax) 1632{ 1633 mydegremin = degreemin; 1634 mydegremax = degreemax; 1635} 1636 1637 1638void Approx_ComputeLine::SetTolerances(const Standard_Real Tolerance3d, 1639 const Standard_Real Tolerance2d) 1640{ 1641 mytol3d = Tolerance3d; 1642 mytol2d = Tolerance2d; 1643} 1644 1645 1646void Approx_ComputeLine::SetConstraints(const AppParCurves_Constraint FirstC, 1647 const AppParCurves_Constraint LastC) 1648{ 1649 myfirstC = FirstC; 1650 mylastC = LastC; 1651} 1652 1653 1654 1655Standard_Boolean Approx_ComputeLine::IsAllApproximated() 1656const { 1657 return alldone; 1658} 1659 1660Standard_Boolean Approx_ComputeLine::IsToleranceReached() 1661const { 1662 return tolreached; 1663} 1664 1665void Approx_ComputeLine::Error(const Standard_Integer Index, 1666 Standard_Real& tol3d, 1667 Standard_Real& tol2d) const 1668{ 1669 tol3d = Tolers3d.Value(Index); 1670 tol2d = Tolers2d.Value(Index); 1671} 1672 1673Approx_ParametrizationType Approx_ComputeLine::Parametrization() const 1674{ 1675 return Par; 1676} 1677