src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DEdge.cxx

   1 // Copyright (C) 2007-2016  CEA/DEN, EDF R&D
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License, or (at your option) any later version.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19 // Author : Anthony Geay (CEA/DEN)
  20
  21 #include "InterpKernelGeo2DEdge.hxx"
  22 #include "InterpKernelGeo2DEdgeLin.hxx"
  23 #include "InterpKernelGeo2DEdgeInfLin.hxx"
  24 //#include "EdgeParabol.hxx"
  25 #include "InterpKernelGeo2DEdgeArcCircle.hxx"
  26 #include "InterpKernelGeo2DComposedEdge.hxx"
  27 #include "InterpKernelException.hxx"
  28
  29 #include <algorithm>
  30
  31 using namespace INTERP_KERNEL;
  32
  33 MergePoints::MergePoints():_ass1Start1(0),_ass1End1(0),_ass1Start2(0),_ass1End2(0),
  34     _ass2Start1(0),_ass2End1(0),_ass2Start2(0),_ass2End2(0)
  35 {
  36 }
  37
  38 void MergePoints::start1Replaced()
  39 {
  40   unsigned nbOfAsso=getNumberOfAssociations();
  41   if(nbOfAsso==0)
  42     _ass1Start1=1;
  43   else
  44     _ass2Start1=1;
  45 }
  46
  47 void MergePoints::end1Replaced()
  48 {
  49   unsigned nbOfAsso=getNumberOfAssociations();
  50   if(nbOfAsso==0)
  51     _ass1End1=1;
  52   else
  53     _ass2End1=1;
  54 }
  55
  56 void MergePoints::start1OnStart2()
  57 {
  58   unsigned nbOfAsso=getNumberOfAssociations();
  59   if(nbOfAsso==0)
  60     {
  61       _ass1Start1=1;
  62       _ass1Start2=1;
  63     }
  64   else
  65     {
  66       _ass2Start1=1;
  67       _ass2Start2=1;
  68     }
  69 }
  70
  71 void MergePoints::start1OnEnd2()
  72 {
  73   unsigned nbOfAsso=getNumberOfAssociations();
  74   if(nbOfAsso==0)
  75     {
  76       _ass1Start1=1;
  77       _ass1End2=1;
  78     }
  79   else
  80     {
  81       _ass2Start1=1;
  82       _ass2End2=1;
  83     }
  84 }
  85
  86 void MergePoints::end1OnStart2()
  87 {
  88   unsigned nbOfAsso=getNumberOfAssociations();
  89   if(nbOfAsso==0)
  90     {
  91       _ass1End1=1;
  92       _ass1Start2=1;
  93     }
  94   else
  95     {
  96       _ass2End1=1;
  97       _ass2Start2=1;
  98     }
  99 }
 100
 101 void MergePoints::end1OnEnd2()
 102 {
 103   unsigned nbOfAsso=getNumberOfAssociations();
 104   if(nbOfAsso==0)
 105     {
 106       _ass1End1=1;
 107       _ass1End2=1;
 108     }
 109   else
 110     {
 111       _ass2End1=1;
 112       _ass2End2=1;
 113     }
 114 }
 115
 116 bool MergePoints::isStart1(unsigned rk) const
 117 {
 118   if(rk==0)
 119     return _ass1Start1;
 120   else
 121     return _ass2Start1;
 122 }
 123
 124 bool MergePoints::isEnd1(unsigned rk) const
 125 {
 126   if(rk==0)
 127     return _ass1End1;
 128   else
 129     return _ass2End1;
 130 }
 131
 132 bool MergePoints::isStart2(unsigned rk) const
 133 {
 134   if(rk==0)
 135     return _ass1Start2;
 136   else
 137     return _ass2Start2;
 138 }
 139
 140 bool MergePoints::isEnd2(unsigned rk) const
 141 {
 142   if(rk==0)
 143     return _ass1End2;
 144   else
 145     return _ass2End2;
 146 }
 147
 148 void MergePoints::clear()
 149 {
 150   _ass1Start1=0;_ass1End1=0;_ass1Start2=0;_ass1End2=0;
 151   _ass2Start1=0;_ass2End1=0;_ass2Start2=0;_ass2End2=0;
 152 }
 153
 154 unsigned MergePoints::getNumberOfAssociations() const
 155 {
 156   unsigned ret=0;
 157   unsigned subTot=_ass1Start1+_ass1End1+_ass1Start2+_ass1End2;
 158   if(subTot!=0)
 159     ret++;
 160   subTot=_ass2Start1+_ass2End1+_ass2Start2+_ass2End2;
 161   if(subTot!=0)
 162     ret++;
 163   return ret;
 164 }
 165
 166 void MergePoints::PushInMap(int key, int value, std::map<int,int>& mergedNodes)
 167 {
 168   if(key!=-1 && value!=-1)
 169     mergedNodes[key]=value;
 170 }
 171
 172 void MergePoints::updateMergedNodes(int e1Start, int e1End, int e2Start, int e2End, std::map<int,int>& mergedNodes)
 173 {
 174   unsigned subTot(_ass1Start1+_ass1End1+_ass1Start2+_ass1End2);
 175   if(subTot!=0)
 176     {
 177       if(_ass1Start1 && _ass1Start2)
 178         PushInMap(e2Start,e1Start,mergedNodes);
 179       if(_ass1Start1 && _ass1End2)
 180         PushInMap(e2End,e1Start,mergedNodes);
 181       if(_ass1End1 && _ass1Start2)
 182         PushInMap(e2Start,e1End,mergedNodes);
 183       if(_ass1End1 && _ass1End2)
 184         PushInMap(e2End,e1End,mergedNodes);
 185     }
 186   subTot=_ass2Start1+_ass2End1+_ass2Start2+_ass2End2;
 187   if(subTot!=0)
 188     {
 189       if(_ass2Start1 && _ass2Start2)
 190         PushInMap(e2Start,e1Start,mergedNodes);
 191       if(_ass2Start1 && _ass2End2)
 192         PushInMap(e2End,e1Start,mergedNodes);
 193       if(_ass2End1 && _ass2Start2)
 194         PushInMap(e2Start,e1End,mergedNodes);
 195       if(_ass2End1 && _ass2End2)
 196         PushInMap(e2End,e1End,mergedNodes);
 197     }
 198 }
 199
 200 IntersectElement::IntersectElement(double val1, double val2, bool start1, bool end1, bool start2, bool end2, Node *node
 201                                    , const Edge& e1, const Edge& e2, bool keepOrder):_1S(keepOrder?start1:start2),
 202                                        _1E(keepOrder?end1:end2),
 203                                        _2S(keepOrder?start2:start1),
 204                                        _2E(keepOrder?end2:end1),
 205                                        _chararct_val_for_e1(keepOrder?val1:val2),
 206                                        _chararct_val_for_e2(keepOrder?val2:val1),
 207                                        _node(node),_loc_of_node(node->getLoc()),_e1(keepOrder?e1:e2),
 208                                        _e2(keepOrder?e2:e1)
 209 {
 210 }
 211
 212 IntersectElement::IntersectElement(const IntersectElement& other):_1S(other._1S),_1E(other._1E),_2S(other._2S),_2E(other._2E),
 213     _chararct_val_for_e1(other._chararct_val_for_e1),
 214     _chararct_val_for_e2(other._chararct_val_for_e2),_node(other._node),
 215     _loc_of_node(other._loc_of_node),_e1(other._e1), _e2(other._e2)
 216 {
 217   if(_node)
 218     _node->incrRef();
 219 }
 220
 221 IntersectElement& IntersectElement::operator=(const IntersectElement& other)
 222 {
 223   _1S=other._1S;_1E=other._1E; _2S=other._2S; _2E=other._2E;
 224   _chararct_val_for_e1=other._chararct_val_for_e1;
 225   _chararct_val_for_e2=other._chararct_val_for_e2;
 226   setNode(other._node);
 227   return *this;
 228 }
 229
 230 bool IntersectElement::operator<(const IntersectElement& other) const
 231 {
 232   return _e1.isLower(_chararct_val_for_e1,other._chararct_val_for_e1);
 233 }
 234
 235 IntersectElement::~IntersectElement()
 236 {
 237   if(_node)
 238     _node->decrRef();
 239 }
 240
 241 /*!
 242  * Returns 0 or 1.
 243  */
 244 bool IntersectElement::isOnMergedExtremity() const
 245 {
 246   if( (_1S && _2S) || (_1S && _2E) || (_1E && _2S) || (_1E && _2E) )
 247     return true;
 248   return false;
 249 }
 250
 251 /*!
 252  * To call if isOnMergedExtremity returned true.
 253  */
 254 void IntersectElement::performMerging(MergePoints& commonNode) const
 255 {
 256   if(_1S && _2S)
 257     {
 258       if(_e1.changeStartNodeWith(_e2.getStartNode()))
 259         {
 260           _e2.getStartNode()->declareOnLim();
 261           commonNode.start1OnStart2();
 262         }
 263     }
 264   else if(_1S && _2E)
 265     {
 266       if(_e1.changeStartNodeWith(_e2.getEndNode()))
 267         {
 268           _e2.getEndNode()->declareOnLim();
 269           commonNode.start1OnEnd2();
 270         }
 271     }
 272   else if(_1E && _2S)
 273     {
 274       if(_e1.changeEndNodeWith(_e2.getStartNode()))
 275         {
 276           _e2.getStartNode()->declareOnLim();
 277           commonNode.end1OnStart2();
 278         }
 279     }
 280   else if(_1E && _2E)
 281     {
 282       if(_e1.changeEndNodeWith(_e2.getEndNode()))
 283         {
 284           _e2.getEndNode()->declareOnLim();
 285           commonNode.end1OnEnd2();
 286         }
 287     }
 288 }
 289
 290 /*!
 291  * This methode is const because 'node' is supposed to be equal geomitrically to _node.
 292  */
 293 void IntersectElement::setNode(Node *node) const
 294 {
 295   if(node!=_node)
 296     {
 297       if(_node)
 298         ((Node *)_node)->decrRef();
 299       (const_cast<IntersectElement *>(this))->_node=node;
 300       if(_node)
 301         _node->incrRef();
 302     }
 303 }
 304
 305 bool IntersectElement::isLowerOnOther(const IntersectElement& other) const
 306 {
 307   return _e2.isLower(_chararct_val_for_e2,other._chararct_val_for_e2);
 308 }
 309
 310 unsigned IntersectElement::isOnExtrForAnEdgeAndInForOtherEdge() const
 311 {
 312   if(( _1S && !(_2S || _2E) ) || ( _1E && !(_2S || _2E) ))
 313     {
 314       if(_1S && !(_2S || _2E))
 315         setNode(_e1.getStartNode());
 316       else
 317         setNode(_e1.getEndNode());
 318       if(_e2.isIn(_chararct_val_for_e2))
 319         return LIMIT_ON;
 320       return LIMIT_ALONE;
 321     }
 322   if(( _2S && !(_1S || _1E) ) || ( _2E && !(_1S || _1E)))
 323     {
 324       if(_2S && !(_1S || _1E))
 325         setNode(_e2.getStartNode());
 326       else
 327         setNode(_e2.getEndNode());
 328       if(_e1.isIn(_chararct_val_for_e1))
 329         return LIMIT_ON;
 330       return LIMIT_ALONE;
 331     }
 332   return NO_LIMIT;
 333 }
 334
 335 bool IntersectElement::isIncludedByBoth() const
 336 {
 337   return _e1.isIn(_chararct_val_for_e1) && _e2.isIn(_chararct_val_for_e2);
 338 }
 339
 340 bool EdgeIntersector::intersect(const Bounds *whereToFind, std::vector<Node *>& newNodes, bool& order, MergePoints& commonNode)
 341 {
 342   std::list< IntersectElement > listOfIntesc=getIntersectionsCharacteristicVal();
 343   std::list< IntersectElement >::iterator iter;
 344   for(iter=listOfIntesc.begin();iter!=listOfIntesc.end();)
 345     {
 346       if((*iter).isOnMergedExtremity())
 347         {
 348           (*iter).performMerging(commonNode);
 349           iter=listOfIntesc.erase(iter);
 350           continue;
 351         }
 352       unsigned tmp=(*iter).isOnExtrForAnEdgeAndInForOtherEdge();
 353       if(tmp==IntersectElement::LIMIT_ALONE)
 354         {
 355           iter=listOfIntesc.erase(iter);
 356           continue;
 357         }
 358       else if(tmp==IntersectElement::LIMIT_ON)
 359         {
 360           (*iter).attachLoc();
 361           iter++;
 362           continue;
 363         }
 364       if(!(*iter).isIncludedByBoth())
 365         {
 366           iter=listOfIntesc.erase(iter);
 367           continue;
 368         }
 369       (*iter).attachLoc();
 370       iter++;
 371     }
 372   if(listOfIntesc.size()==0)
 373     return false;
 374   if(listOfIntesc.size()==1)
 375     {
 376       order=true;//useless
 377       newNodes.push_back(listOfIntesc.front().getNodeAndReleaseIt());
 378     }
 379   else
 380     {
 381       std::vector<IntersectElement> vecOfIntesc(listOfIntesc.begin(),listOfIntesc.end());
 382       listOfIntesc.clear();
 383       sort(vecOfIntesc.begin(),vecOfIntesc.end());
 384       for(std::vector<IntersectElement>::iterator iterV=vecOfIntesc.begin();iterV!=vecOfIntesc.end();iterV++)
 385         newNodes.push_back((*iterV).getNodeAndReleaseIt());
 386       order=vecOfIntesc.front().isLowerOnOther(vecOfIntesc.back());
 387     }
 388   return true;
 389 }
 390
 391 /*!
 392  * Locates 'node' regarding edge this->_e1. If node is located close to (with distant lt epsilon) start or end point of _e1,
 393  * 'node' takes its place. In this case 'obvious' is set to true and 'commonNode' stores information of merge point and finally 'where' is set.
 394  * Furthermore 'node' is declared as ON LIMIT to indicate in locating process that an absolute location computation will have to be done.
 395  * If 'node' is not close to start or end point of _e1, 'obvious' is set to false and 'commonNode' and 'where' are let unchanged.
 396  */
 397 void EdgeIntersector::obviousCaseForCurvAbscisse(Node *node, TypeOfLocInEdge& where, MergePoints& commonNode, bool& obvious) const
 398 {
 399   obvious=true;
 400   if(node->isEqual(*_e1.getStartNode()))
 401     {
 402       where=START;
 403       if(_e1.changeStartNodeWith(node))
 404         {
 405           commonNode.start1Replaced();
 406           node->declareOnLim();
 407         }
 408       return ;
 409     }
 410   if(node->isEqual(*_e1.getEndNode()))
 411     {
 412       where=END;
 413       if(_e1.changeEndNodeWith(node))
 414         {
 415           commonNode.end1Replaced();
 416           node->declareOnLim();
 417         }
 418       return ;
 419     }
 420   obvious=false;
 421 }
 422
 423 Edge::Edge(double sX, double sY, double eX, double eY):_cnt(1),_loc(FULL_UNKNOWN),_start(new Node(sX,sY)),_end(new Node(eX,eY))
 424 {
 425 }
 426
 427 Edge::~Edge()
 428 {
 429   _start->decrRef();
 430   if(_end)
 431     _end->decrRef();
 432 }
 433
 434 bool Edge::decrRef()
 435 {
 436   bool ret=(--_cnt==0);
 437   if(ret)
 438     delete this;
 439   return ret;
 440 }
 441
 442 void Edge::declareOn() const
 443 {
 444   if(_loc==FULL_UNKNOWN)
 445     {
 446       _loc=FULL_ON_1;
 447       _start->declareOn();
 448       _end->declareOn();
 449     }
 450 }
 451
 452 void Edge::declareIn() const
 453 {
 454   if(_loc==FULL_UNKNOWN)
 455     {
 456       _loc=FULL_IN_1;
 457       _start->declareIn();
 458       _end->declareIn();
 459     }
 460 }
 461
 462 void Edge::declareOut() const
 463 {
 464   if(_loc==FULL_UNKNOWN)
 465     {
 466       _loc=FULL_OUT_1;
 467       _start->declareOut();
 468       _end->declareOut();
 469     }
 470 }
 471
 472 void Edge::fillXfigStreamForLoc(std::ostream& stream) const
 473 {
 474   switch(_loc)
 475   {
 476     case FULL_IN_1:
 477       stream << '2';//Green
 478       break;
 479     case FULL_OUT_1:
 480       stream << '1';//Bleue
 481       break;
 482     case FULL_ON_1:
 483       stream << '4';//Red
 484       break;
 485     default:
 486       stream << '0';
 487   }
 488 }
 489
 490 bool Edge::changeStartNodeWith(Node *otherStartNode) const
 491 {
 492   if(_start==otherStartNode)
 493     return true;
 494   if(_start->isEqual(*otherStartNode))
 495     {
 496       ((const_cast<Edge *>(this))->_start)->decrRef();//un-const cast Ok thanks to 2 lines above.
 497       ((const_cast<Edge *>(this))->_start)=otherStartNode;
 498       _start->incrRef();
 499       return true;
 500     }
 501   return false;
 502 }
 503
 504 bool Edge::changeStartNodeWithAndKeepTrack(Node *otherStartNode, std::vector<Node *>& track) const
 505 {
 506   if(_start==otherStartNode)
 507     return true;
 508   if(_start->isEqualAndKeepTrack(*otherStartNode,track))
 509     {
 510       ((const_cast<Edge *>(this))->_start)->decrRef();//un-const cast Ok thanks to 2 lines above.
 511       ((const_cast<Edge *>(this))->_start)=otherStartNode;
 512       otherStartNode->incrRef();
 513       return true;
 514     }
 515   return false;
 516 }
 517
 518 bool Edge::changeEndNodeWith(Node *otherEndNode) const
 519 {
 520   if(_end==otherEndNode)
 521     return true;
 522   if(_end->isEqual(*otherEndNode))
 523     {
 524       ((const_cast<Edge *>(this))->_end)->decrRef();
 525       ((const_cast<Edge *>(this))->_end)=otherEndNode;
 526       _end->incrRef();
 527       return true;
 528     }
 529   return false;
 530 }
 531
 532 bool Edge::changeEndNodeWithAndKeepTrack(Node *otherEndNode, std::vector<Node *>& track) const
 533 {
 534   if(_end==otherEndNode)
 535     return true;
 536   if(_end->isEqualAndKeepTrack(*otherEndNode,track))
 537     {
 538       ((const_cast<Edge *>(this))->_end)->decrRef();
 539       ((const_cast<Edge *>(this))->_end)=otherEndNode;
 540       otherEndNode->incrRef();
 541       return true;
 542     }
 543   return false;
 544 }
 545
 546 /*!
 547  * Precondition : 'start' and 'end' are lying on the same curve than 'this'.
 548  * Add in vec the sub edge lying on this.
 549  * If 'start' is equal (by pointer) to '_end' and 'end' is equal to '_end' too nothing is added.
 550  * If 'start' is equal (by pointer) to '_start' and 'end' is equal to '_start' too nothing is added.
 551  * If 'start' is equal (by pointer) to '_start' and 'end' is equal to '_end' this is added in vec.
 552  */
 553 void Edge::addSubEdgeInVector(Node *start, Node *end, ComposedEdge& vec) const
 554 {
 555   if((start==_start && end==_start) || (start==_end && end==_end))
 556     return ;
 557   if(start==_start && end==_end)
 558     {
 559       incrRef();
 560       vec.pushBack(const_cast<Edge *>(this));
 561       return ;
 562     }
 563   vec.pushBack(buildEdgeLyingOnMe(start,end,true));
 564 }
 565
 566 /*!
 567  * Retrieves a vector 'vectOutput' that is normal to 'this'. 'vectOutput' is normalized.
 568  */
 569 void Edge::getNormalVector(double *vectOutput) const
 570 {
 571   std::copy((const double *)(*_end),(const double *)(*_end)+2,vectOutput);
 572   std::transform(vectOutput,vectOutput+2,(const double *)(*_start),vectOutput,std::minus<double>());
 573   double norm=1./Node::norm(vectOutput);
 574   std::transform(vectOutput,vectOutput+2,vectOutput,bind2nd(std::multiplies<double>(),norm));
 575   double tmp=vectOutput[0];
 576   vectOutput[0]=vectOutput[1];
 577   vectOutput[1]=-tmp;
 578 }
 579
 580 Edge *Edge::BuildEdgeFrom3Points(const double *start, const double *middle, const double *end)
 581 {
 582   Node *b(new Node(start[0],start[1])),*m(new Node(middle[0],middle[1])),*e(new Node(end[0],end[1]));
 583   EdgeLin *e1(new EdgeLin(b,m)),*e2(new EdgeLin(m,e));
 584   SegSegIntersector inters(*e1,*e2); bool colinearity=inters.areColinears(); delete e1; delete e2;
 585   Edge *ret=0;
 586   if(colinearity)
 587     ret=new EdgeLin(b,e);
 588   else
 589     ret=new EdgeArcCircle(b,m,e);
 590   b->decrRef(); m->decrRef(); e->decrRef();
 591   return ret;
 592 }
 593
 594 Edge *Edge::BuildEdgeFrom(Node *start, Node *end)
 595 {
 596   return new EdgeLin(start,end);
 597 }
 598
 599 Edge *Edge::BuildFromXfigLine(std::istream& str)
 600 {
 601   unsigned char type;
 602   str >> type;
 603   if(type=='2')
 604     return new EdgeLin(str);
 605   else if(type=='5')
 606     return new EdgeArcCircle(str);
 607   else
 608     {
 609       std::cerr << "Unknown line found...";
 610       return 0;
 611     }
 612 }
 613
 614 /*!
 615  * \param other The Edge with which we are going to intersect.
 616  * \param commonNode Output. The common nodes found during operation of intersecting.
 617  * \param outVal1 Output filled in case true is returned. It specifies the new or not new edges by which 'this' is replaced after intersecting op.
 618  * \param outVal2 Output filled in case true is returned. It specifies the new or not new edges by which 'other' is replaced after intersecting op.
 619  * return true if the intersection between this.
 620  */
 621 bool Edge::intersectWith(const Edge *other, MergePoints& commonNode,
 622                          ComposedEdge& outVal1, ComposedEdge& outVal2) const
 623 {
 624   bool ret=true;
 625   Bounds *merge=_bounds.nearlyAmIIntersectingWith(other->getBounds());
 626   if(!merge)
 627     return false;
 628   delete merge;
 629   merge=0;
 630   EdgeIntersector *intersector=BuildIntersectorWith(this,other);
 631   ret=Intersect(this,other,intersector,merge,commonNode,outVal1,outVal2);
 632   delete intersector;
 633   return ret;
 634 }
 635
 636 bool Edge::IntersectOverlapped(const Edge *f1, const Edge *f2, EdgeIntersector *intersector, MergePoints& commonNode,
 637                                ComposedEdge& outValForF1, ComposedEdge& outValForF2)
 638 {
 639   bool rev=intersector->haveTheySameDirection();
 640   Node *f2Start=f2->getNode(rev?START:END);
 641   Node *f2End=f2->getNode(rev?END:START);
 642   TypeOfLocInEdge place1, place2;
 643   intersector->getPlacements(f2Start,f2End,place1,place2,commonNode);
 644   int codeForIntersectionCase=CombineCodes(place1,place2);
 645   return SplitOverlappedEdges(f1,f2,f2Start,f2End,rev,codeForIntersectionCase,outValForF1,outValForF2);
 646 }
 647
 648 /*!
 649  * Perform 1D linear interpolation. Warning distrib1 and distrib2 are expected to be in ascending mode.
 650  */
 651 void Edge::Interpolate1DLin(const std::vector<double>& distrib1, const std::vector<double>& distrib2, std::map<int, std::map<int,double> >& result)
 652 {
 653   int nbOfV1=distrib1.size()-1;
 654   int nbOfV2=distrib2.size()-1;
 655   Node *n1=new Node(0.,0.); Node *n3=new Node(0.,0.);
 656   Node *n2=new Node(0.,0.); Node *n4=new Node(0.,0.);
 657   MergePoints commonNode;
 658   for(int i=0;i<nbOfV1;i++)
 659     {
 660       std::vector<double>::const_iterator iter=find_if(distrib2.begin()+1,distrib2.end(),bind2nd(std::greater_equal<double>(),distrib1[i]));
 661       if(iter!=distrib2.end())
 662         {
 663           for(int j=(iter-1)-distrib2.begin();j<nbOfV2;j++)
 664             {
 665               if(distrib2[j]<=distrib1[i+1])
 666                 {
 667                   EdgeLin *e1=new EdgeLin(n1,n2); EdgeLin *e2=new EdgeLin(n3,n4);
 668                   n1->setNewCoords(distrib1[i],0.); n2->setNewCoords(distrib1[i+1],0.);
 669                   n3->setNewCoords(distrib2[j],0.); n4->setNewCoords(distrib2[j+1],0.);
 670                   ComposedEdge *f1=new ComposedEdge;
 671                   ComposedEdge *f2=new ComposedEdge;
 672                   SegSegIntersector inters(*e1,*e2);
 673                   bool b1,b2;
 674                   inters.areOverlappedOrOnlyColinears(0,b1,b2);
 675                   if(IntersectOverlapped(e1,e2,&inters,commonNode,*f1,*f2))
 676                     {
 677                       result[i][j]=f1->getCommonLengthWith(*f2)/e1->getCurveLength();
 678                     }
 679                   ComposedEdge::Delete(f1); ComposedEdge::Delete(f2);
 680                   e1->decrRef(); e2->decrRef();
 681                 }
 682             }
 683         }
 684     }
 685   n1->decrRef(); n2->decrRef(); n3->decrRef(); n4->decrRef();
 686 }
 687
 688 EdgeIntersector *Edge::BuildIntersectorWith(const Edge *e1, const Edge *e2)
 689 {
 690   EdgeIntersector *ret=0;
 691   const EdgeLin *tmp1=0;
 692   const EdgeArcCircle *tmp2=0;
 693   unsigned char type1=e1->getTypeOfFunc();
 694   e1->dynCastFunction(tmp1,tmp2);
 695   unsigned char type2=e2->getTypeOfFunc();
 696   e2->dynCastFunction(tmp1,tmp2);
 697   type1|=type2;
 698   switch(type1)
 699   {
 700     case 1:// Intersection seg/seg
 701       ret=new SegSegIntersector((const EdgeLin &)(*e1),(const EdgeLin &)(*e2));
 702       break;
 703     case 5:// Intersection seg/arc of circle
 704       ret=new ArcCSegIntersector(*tmp2,*tmp1,tmp2==e1);
 705       break;
 706     case 4:// Intersection arc/arc of circle
 707       ret=new ArcCArcCIntersector((const EdgeArcCircle &)(*e1),(const EdgeArcCircle &)(*e2));
 708       break;
 709     default:
 710       //Should never happen
 711       throw Exception("A non managed association of edge has been detected. Go work for intersection computation implementation.");
 712   }
 713   return ret;
 714 }
 715
 716 /*!
 717  * See Node::applySimilarity to see signification of params.
 718  */
 719 void Edge::applySimilarity(double xBary, double yBary, double dimChar)
 720 {
 721   _bounds.applySimilarity(xBary,yBary,dimChar);
 722 }
 723
 724 void Edge::unApplySimilarity(double xBary, double yBary, double dimChar)
 725 {
 726   _bounds.unApplySimilarity(xBary,yBary,dimChar);
 727 }
 728
 729 void Edge::getMiddleOfPointsOriented(const double *p1, const double *p2, double *mid) const
 730 {
 731   return getMiddleOfPoints(p1, p2, mid);
 732 }
 733
 734 bool Edge::Intersect(const Edge *f1, const Edge *f2, EdgeIntersector *intersector, const Bounds *whereToFind, MergePoints& commonNode,
 735                      ComposedEdge& outValForF1, ComposedEdge& outValForF2)
 736 {
 737   bool obviousNoIntersection;
 738   bool areOverlapped;
 739   intersector->areOverlappedOrOnlyColinears(whereToFind,obviousNoIntersection,areOverlapped);
 740   if(areOverlapped)
 741     return IntersectOverlapped(f1,f2,intersector,commonNode,outValForF1,outValForF2);
 742   if(obviousNoIntersection)
 743     return false;
 744   std::vector<Node *> newNodes;
 745   bool order;
 746   if(intersector->intersect(whereToFind,newNodes,order,commonNode))
 747     {
 748       if(newNodes.empty())
 749         throw Exception("Internal error occured - error in intersector implementation!");// This case should never happen
 750       std::vector<Node *>::iterator iter=newNodes.begin();
 751       std::vector<Node *>::reverse_iterator iterR=newNodes.rbegin();
 752       f1->addSubEdgeInVector(f1->getStartNode(),*iter,outValForF1);
 753       f2->addSubEdgeInVector(f2->getStartNode(),order?*iter:*iterR,outValForF2);
 754       for(std::vector<Node *>::iterator iter2=newNodes.begin();iter2!=newNodes.end();iter2++,iterR++)
 755         {
 756           if((iter2+1)==newNodes.end())
 757             {
 758               f1->addSubEdgeInVector(*iter2,f1->getEndNode(),outValForF1);
 759               (*iter2)->decrRef();
 760               f2->addSubEdgeInVector(order?*iter2:*iterR,f2->getEndNode(),outValForF2);
 761             }
 762           else
 763             {
 764               f1->addSubEdgeInVector(*iter2,*(iter2+1),outValForF1);
 765               (*iter2)->decrRef();
 766               f2->addSubEdgeInVector(order?*iter2:*iterR,order?*(iter2+1):*(iterR+1),outValForF2);
 767             }
 768         }
 769       return true;
 770     }
 771   else//no intersection inside whereToFind
 772     return false;
 773 }
 774
 775 int Edge::CombineCodes(TypeOfLocInEdge code1, TypeOfLocInEdge code2)
 776 {
 777   int ret=(int)code1;
 778   ret*=OFFSET_FOR_TYPEOFLOCINEDGE;
 779   ret+=(int)code2;
 780   return ret;
 781 }
 782
 783 /*!
 784  * This method splits e1 and e2 into pieces as much sharable as possible. The precondition to the call of this method
 785  * is that e1 and e2 have been declared as overlapped by corresponding intersector built from e1 and e2 type.
 786  *
 787  * @param nS start node of e2 with the SAME DIRECTION as e1. The pointer nS should be equal to start node of e2 or to its end node.
 788  * @param nE end node of e2 with the SAME DIRECTION as e1. The pointer nE should be equal to start node of e2 or to its end node.
 789  * @param direction is param that specifies if e2 and e1 have same directions (true) or opposed (false).
 790  * @param code is the code returned by method Edge::combineCodes.
 791  */
 792 bool Edge::SplitOverlappedEdges(const Edge *e1, const Edge *e2, Node *nS, Node *nE, bool direction, int code,
 793                                 ComposedEdge& outVal1, ComposedEdge& outVal2)
 794 {
 795   Edge *tmp;
 796   switch(code)
 797   {
 798     case OUT_BEFORE*OFFSET_FOR_TYPEOFLOCINEDGE+START:      // OUT_BEFORE - START
 799     case OUT_BEFORE*OFFSET_FOR_TYPEOFLOCINEDGE+OUT_BEFORE: // OUT_BEFORE - OUT_BEFORE
 800     case OUT_AFTER*OFFSET_FOR_TYPEOFLOCINEDGE+OUT_AFTER:   // OUT_AFTER - OUT_AFTER
 801     case END*OFFSET_FOR_TYPEOFLOCINEDGE+OUT_AFTER:         // END - OUT_AFTER
 802     case END*OFFSET_FOR_TYPEOFLOCINEDGE+START:             // END - START
 803     return false;
 804     case INSIDE*OFFSET_FOR_TYPEOFLOCINEDGE+OUT_AFTER:      // INSIDE - OUT_AFTER
 805     outVal1.pushBack(e1->buildEdgeLyingOnMe(e1->getStartNode(),nS,true));
 806     tmp=e1->buildEdgeLyingOnMe(nS,e1->getEndNode()); tmp->incrRef();
 807     outVal1.pushBack(tmp);
 808     outVal2.resize(2);
 809     outVal2.setValueAt(direction?0:1,tmp,direction); tmp->declareOn();
 810     outVal2.setValueAt(direction?1:0,e1->buildEdgeLyingOnMe(e1->getEndNode(),nE,direction));
 811     return true;
 812     case INSIDE*OFFSET_FOR_TYPEOFLOCINEDGE+INSIDE:         // INSIDE - INSIDE
 813     {
 814       if(!e2->isIn(e2->getCharactValue(*(e1->getStartNode()))))
 815         {
 816           e2->incrRef(); e2->incrRef();
 817           outVal1.resize(3);
 818           outVal1.setValueAt(0,e1->buildEdgeLyingOnMe(e1->getStartNode(),nS));
 819           outVal1.setValueAt(1,const_cast<Edge*>(e2),direction);
 820           outVal1.setValueAt(2,e1->buildEdgeLyingOnMe(nE,e1->getEndNode()));
 821           outVal2.pushBack(const_cast<Edge*>(e2)); e2->declareOn();
 822           return true;
 823         }
 824       else
 825         {
 826           outVal1.resize(3);
 827           outVal2.resize(3);
 828           tmp=e1->buildEdgeLyingOnMe(e1->getStartNode(),nE); tmp->incrRef(); tmp->declareOn();
 829           outVal1.setValueAt(0,tmp,true); outVal2.setValueAt(direction?2:0,tmp,direction);
 830           outVal1.setValueAt(1,e1->buildEdgeLyingOnMe(nE,nS));
 831           tmp=e1->buildEdgeLyingOnMe(nS,e1->getEndNode()); tmp->incrRef(); tmp->declareOn();
 832           outVal1.setValueAt(2,tmp,true); outVal2.setValueAt(direction?0:2,tmp,direction);
 833           tmp=e1->buildEdgeLyingOnMe(e1->getEndNode(),e1->getStartNode());
 834           outVal2.setValueAt(1,tmp,direction);
 835           return true;
 836         }
 837     }
 838     case OUT_BEFORE*OFFSET_FOR_TYPEOFLOCINEDGE+INSIDE:     // OUT_BEFORE - INSIDE
 839     {
 840       tmp=e1->buildEdgeLyingOnMe(e1->getStartNode(),nE); tmp->incrRef();
 841       outVal1.pushBack(tmp);
 842       outVal1.pushBack(e1->buildEdgeLyingOnMe(nE,e1->getEndNode()));
 843       outVal2.resize(2);
 844       outVal2.setValueAt(direction?0:1,e1->buildEdgeLyingOnMe(nS,e1->getStartNode(),direction));
 845       outVal2.setValueAt(direction?1:0,tmp,direction); tmp->declareOn();
 846       return true;
 847     }
 848     case OUT_BEFORE*OFFSET_FOR_TYPEOFLOCINEDGE+OUT_AFTER:  // OUT_BEFORE - OUT_AFTER
 849     {
 850       e1->incrRef(); e1->incrRef();
 851       outVal1.pushBack(const_cast<Edge*>(e1));
 852       outVal2.resize(3);
 853       outVal2.setValueAt(direction?0:2,e1->buildEdgeLyingOnMe(nS,e1->getStartNode(),direction));
 854       outVal2.setValueAt(1,const_cast<Edge*>(e1),direction); e1->declareOn();
 855       outVal2.setValueAt(direction?2:0,e1->buildEdgeLyingOnMe(e1->getEndNode(),nE,direction));
 856       return true;
 857     }
 858     case START*OFFSET_FOR_TYPEOFLOCINEDGE+END:             // START - END
 859     {
 860       e1->incrRef(); e1->incrRef();
 861       outVal1.pushBack(const_cast<Edge*>(e1));
 862       outVal2.pushBack(const_cast<Edge*>(e1),direction); e1->declareOn();
 863       return true;
 864     }
 865     case START*OFFSET_FOR_TYPEOFLOCINEDGE+OUT_AFTER:       // START - OUT_AFTER
 866     {
 867       e1->incrRef(); e1->incrRef();
 868       outVal1.pushBack(const_cast<Edge*>(e1));
 869       outVal2.resize(2);
 870       outVal2.setValueAt(direction?0:1,const_cast<Edge*>(e1),direction); e1->declareOn();
 871       outVal2.setValueAt(direction?1:0,e1->buildEdgeLyingOnMe(e1->getEndNode(),nE,direction));
 872       return true;
 873     }
 874     case INSIDE*OFFSET_FOR_TYPEOFLOCINEDGE+END:            // INSIDE - END
 875     {
 876       e2->incrRef(); e2->incrRef();
 877       outVal1.pushBack(e1->buildEdgeLyingOnMe(e1->getStartNode(),nS,true));
 878       outVal1.pushBack(const_cast<Edge*>(e2),direction);
 879       outVal2.pushBack(const_cast<Edge*>(e2)); e2->declareOn();
 880       return true;
 881     }
 882     case OUT_BEFORE*OFFSET_FOR_TYPEOFLOCINEDGE+END:        // OUT_BEFORE - END
 883     {
 884       e1->incrRef(); e1->incrRef();
 885       outVal1.pushBack(const_cast<Edge*>(e1));
 886       outVal2.resize(2);
 887       outVal2.setValueAt(direction?0:1,e1->buildEdgeLyingOnMe(nS,e1->getStartNode(),direction));
 888       outVal2.setValueAt(direction?1:0,const_cast<Edge*>(e1),direction); e1->declareOn();
 889       return true;
 890     }
 891     case START*OFFSET_FOR_TYPEOFLOCINEDGE+INSIDE:          // START - INSIDE
 892     {
 893       e2->incrRef(); e2->incrRef();
 894       outVal1.pushBack(const_cast<Edge*>(e2),direction);
 895       outVal1.pushBack(e1->buildEdgeLyingOnMe(nE,e1->getEndNode()));
 896       outVal2.pushBack(const_cast<Edge*>(e2)); e2->declareOn();
 897       return true;
 898     }
 899     case INSIDE*OFFSET_FOR_TYPEOFLOCINEDGE+START:          // INSIDE - START
 900     {
 901       outVal1.resize(2);
 902       outVal2.resize(2);
 903       tmp=e1->buildEdgeLyingOnMe(nS,e1->getEndNode()); tmp->incrRef(); tmp->declareOn();
 904       outVal1.setValueAt(0,e1->buildEdgeLyingOnMe(e1->getStartNode(),nS));
 905       outVal1.setValueAt(1,tmp);
 906       outVal2.setValueAt(direction?0:1,tmp,direction);
 907       outVal2.setValueAt(direction?1:0,e1->buildEdgeLyingOnMe(e1->getEndNode(),nE,direction));
 908       return true;
 909     }
 910     case END*OFFSET_FOR_TYPEOFLOCINEDGE+INSIDE:            // END - INSIDE
 911     {
 912       outVal1.resize(2);
 913       outVal2.resize(2);
 914       tmp=e1->buildEdgeLyingOnMe(e1->getStartNode(),nE); tmp->incrRef(); tmp->declareOn();
 915       outVal1.setValueAt(0,tmp);
 916       outVal1.setValueAt(1,e1->buildEdgeLyingOnMe(nE,e1->getEndNode()));
 917       outVal2.setValueAt(direction?0:1,e1->buildEdgeLyingOnMe(e1->getEndNode(),e1->getStartNode(),direction));
 918       outVal2.setValueAt(direction?1:0,tmp,direction);
 919       return true;
 920     }
 921     default:
 922       throw Exception("Unexpected situation of overlapping edges : internal error occurs ! ");
 923   }
 924 }
 925
 926 bool Edge::isEqual(const Edge& other) const
 927 {
 928   return _start->isEqual(*other._start) && _end->isEqual(*other._end);
 929 }
 930
 931 inline bool eqpair(const std::pair<double,Node *>& p1, const std::pair<double,Node *>& p2)
 932 {
 933   return fabs(p1.first-p2.first)<QUADRATIC_PLANAR::_precision;
 934 }
 935
 936 /**
 937  * This method takes in input nodes in \a subNodes (using \a coo)
 938  *
 939  * \param [in,out] subNodes to be sorted
 940  * \return true if a reordering was necessary false if not.
 941  */
 942 bool Edge::sortSubNodesAbs(const double *coo, std::vector<int>& subNodes)
 943 {
 944   Bounds b;
 945   b.prepareForAggregation();
 946   b.aggregate(getBounds());
 947   double xBary,yBary;
 948   double dimChar(b.getCaracteristicDim());
 949   b.getBarycenter(xBary,yBary);
 950   applySimilarity(xBary,yBary,dimChar);
 951   _start->applySimilarity(xBary,yBary,dimChar);
 952   _end->applySimilarity(xBary,yBary,dimChar);
 953   //
 954   std::size_t sz(subNodes.size()),i(0);
 955   std::vector< std::pair<double,Node *> > an2(sz);
 956   std::map<Node *, int> m;
 957   for(std::vector<int>::const_iterator it=subNodes.begin();it!=subNodes.end();it++,i++)
 958     {
 959       Node *n(new Node(coo[2*(*it)],coo[2*(*it)+1]));
 960       n->applySimilarity(xBary,yBary,dimChar);
 961       m[n]=*it;
 962       an2[i]=std::pair<double,Node *>(getCharactValueBtw0And1(*n),n);
 963     }
 964   std::sort(an2.begin(),an2.end());
 965   //
 966   bool ret(false);
 967   for(i=0;i<sz;i++)
 968     {
 969       int id(m[an2[i].second]);
 970       if(id!=subNodes[i])
 971         { subNodes[i]=id; ret=true; }
 972     }
 973   //
 974   for(std::map<INTERP_KERNEL::Node *,int>::const_iterator it2=m.begin();it2!=m.end();it2++)
 975     (*it2).first->decrRef();
 976   return ret;
 977 }
 978
 979 /**
 980  * Sort nodes so that they all lie consecutively on the edge that has been cut.
 981  */
 982 void Edge::sortIdsAbs(const std::vector<INTERP_KERNEL::Node *>& addNodes, const std::map<INTERP_KERNEL::Node *, int>& mapp1,
 983                       const std::map<INTERP_KERNEL::Node *, int>& mapp2, std::vector<int>& edgesThis)
 984 {
 985   Bounds b;
 986   b.prepareForAggregation();
 987   b.aggregate(getBounds());
 988   double xBary,yBary;
 989   double dimChar=b.getCaracteristicDim();
 990   b.getBarycenter(xBary,yBary);
 991   for(std::vector<Node *>::const_iterator iter=addNodes.begin();iter!=addNodes.end();iter++)
 992     (*iter)->applySimilarity(xBary,yBary,dimChar);
 993   applySimilarity(xBary,yBary,dimChar);
 994   _start->applySimilarity(xBary,yBary,dimChar);
 995   _end->applySimilarity(xBary,yBary,dimChar);
 996   std::size_t sz=addNodes.size();
 997   std::vector< std::pair<double,Node *> > an2(sz);
 998   for(std::size_t i=0;i<sz;i++)
 999     an2[i]=std::pair<double,Node *>(getCharactValueBtw0And1(*addNodes[i]),addNodes[i]);
1000   std::sort(an2.begin(),an2.end());
1001   int startId=(*mapp1.find(_start)).second;
1002   int endId=(*mapp1.find(_end)).second;
1003   std::vector<int> tmpp;
1004   std::vector< std::pair<double,Node *> >::const_iterator itend=std::unique(an2.begin(),an2.end(),eqpair);
1005   for(std::vector< std::pair<double,Node *> >::const_iterator it=an2.begin();it!=itend;it++)
1006     {
1007       int idd=(*mapp2.find((*it).second)).second;
1008       if((*it).first<QUADRATIC_PLANAR::_precision)
1009         {
1010           startId=idd;
1011           continue;
1012         }
1013       if((*it).first>1-QUADRATIC_PLANAR::_precision)
1014         {
1015           endId=idd;
1016           continue;
1017         }
1018       tmpp.push_back(idd);
1019     }
1020   std::vector<int> tmpp2(tmpp.size()+2);
1021   tmpp2[0]=startId;
1022   std::copy(tmpp.begin(),tmpp.end(),tmpp2.begin()+1);
1023   tmpp2[tmpp.size()+1]=endId;
1024   std::vector<int>::iterator itt=std::unique(tmpp2.begin(),tmpp2.end());
1025   tmpp2.resize(std::distance(tmpp2.begin(),itt));
1026   int nbOfEdges=tmpp2.size()-1;
1027   for(int i=0;i<nbOfEdges;i++)
1028     {
1029       edgesThis.push_back(tmpp2[i]);
1030       edgesThis.push_back(tmpp2[i+1]);
1031     }
1032 }