src/HYDROData/HYDROData_Tool.cxx

   1 // Copyright (C) 2014-2015  EDF-R&D
   2 // This library is free software; you can redistribute it and/or
   3 // modify it under the terms of the GNU Lesser General Public
   4 // License as published by the Free Software Foundation; either
   5 // version 2.1 of the License, or (at your option) any later version.
   6 //
   7 // This library is distributed in the hope that it will be useful,
   8 // but WITHOUT ANY WARRANTY; without even the implied warranty of
   9 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  10 // Lesser General Public License for more details.
  11 //
  12 // You should have received a copy of the GNU Lesser General Public
  13 // License along with this library; if not, write to the Free Software
  14 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  15 //
  16 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  17 //
  18
  19 #include <HYDROData_Tool.h>
  20 #include <HYDROData_ArtificialObject.h>
  21 #include <HYDROData_Document.h>
  22 #include <HYDROData_Entity.h>
  23 #include <HYDROData_Iterator.h>
  24 #include <HYDROData_NaturalObject.h>
  25 #include <HYDROData_ShapesGroup.h>
  26 #include <HYDROData_PolylineXY.h>
  27 #include <QColor>
  28 #include <QFile>
  29 #include <QStringList>
  30 #include <QTextStream>
  31 #include <QSet>
  32 #include <BRep_Tool.hxx>
  33 #include <BRepAdaptor_Surface.hxx>
  34 #include <BRepTopAdaptor_FClass2d.hxx>
  35 #include <ElSLib.hxx>
  36 #include <Geom_Curve.hxx>
  37 #include <gp_Pln.hxx>
  38 #include <Quantity_Color.hxx>
  39 #include <TopExp_Explorer.hxx>
  40 #include <TopoDS.hxx>
  41 #include <TopoDS_Face.hxx>
  42 #include <TopoDS_Shape.hxx>
  43 #include <TopoDS_Wire.hxx>
  44 #include <limits>
  45 #include <math.h>
  46
  47
  48 #include <BRepTools.hxx>
  49 #include <NCollection_Map.hxx>
  50 #include <TopTools_ShapeMapHasher.hxx>
  51 #include <BRep_Builder.hxx>
  52 #include <TopTools_IndexedDataMapOfShapeListOfShape.hxx>
  53 #include <TopExp.hxx>
  54 #include <NCollection_List.hxx>
  55 #include <TopTools_ListIteratorOfListOfShape.hxx>
  56
  57 #include <BRepBuilderAPI_MakeFace.hxx>
  58
  59 #include <TopExp_Explorer.hxx>
  60 #include <TopTools_ListIteratorOfListOfShape.hxx>
  61 #include <TopTools_HSequenceOfShape.hxx>
  62
  63 #include <BRep_Builder.hxx>
  64 #include <BRepAlgo_FaceRestrictor.hxx>
  65 #include <BRepCheck_Analyzer.hxx>
  66
  67 #include <ShapeAnalysis.hxx>
  68 #include <ShapeAnalysis_FreeBounds.hxx>
  69
  70 #include <HYDROData_PolylineXY.h>
  71 #include <HYDROData_Polyline3D.h>
  72 #include <HYDROData_Bathymetry.h>
  73 #include <Geom_Line.hxx>
  74 #include <Geom2d_BSplineCurve.hxx>
  75 #include <Standard_Type.hxx>
  76 #include <Geom2d_TrimmedCurve.hxx>
  77 #include <Geom2d_Line.hxx>
  78 #include <Geom_BSplineCurve.hxx>
  79
  80 #include <GeomAPI.hxx>
  81 #include <gp.hxx>
  82 #include <Geom_Plane.hxx>
  83 #include <QFileInfo>
  84 #define BLOCK_SIZE 10000
  85
  86 HYDROData_Tool::ExecStatus HYDROData_Tool::myTriangulationStatus = ExecStatus::None;
  87
  88 static int aMaxNameId = INT_MAX;
  89 static int aMaxColorNb = 92000;
  90 void HYDROData_Tool::WriteStringsToFile( QFile&             theFile,
  91                                          const QStringList& theStrings,
  92                                          const QString&     theSep )
  93 {
  94   if ( !theFile.isOpen() || theStrings.isEmpty() )
  95     return;
  96
  97   QString aWriteStr = theStrings.join( theSep );
  98   if ( aWriteStr.isEmpty() )
  99     return;
 100
 101   QTextStream anOutStream( &theFile );
 102   anOutStream << aWriteStr.toUtf8() << theSep << theSep;
 103 }
 104
 105 bool HYDROData_Tool::ExtractGeneratedObjectName(const QString& theName, int& outValue, QString& thePrefName)
 106 {
 107   QStringList aLs = theName.split('_');
 108   bool ok = false;
 109   QString last = aLs.last();
 110   outValue = last.toInt(&ok);
 111   if (!ok)
 112     return false;
 113   int last_len = last.length();
 114   int total_len = theName.length();
 115   thePrefName = theName.left(total_len-last_len-1);
 116   return true;
 117 }
 118
 119 QString HYDROData_Tool::GenerateObjectName( const Handle(HYDROData_Document)& theDoc,
 120                                             const QString&                    thePrefix,
 121                                             const QStringList&                theUsedNames,
 122                                             const bool                        theIsTryToUsePurePrefix )
 123 {
 124   QStringList aNamesList( theUsedNames );
 125
 126   // Collect all used names in the document
 127   HYDROData_Iterator anIter( theDoc );
 128   for( ; anIter.More(); anIter.Next() )
 129   {
 130     Handle(HYDROData_Entity) anObject = anIter.Current();
 131     if( anObject.IsNull() )
 132       continue;
 133
 134     QString anObjName = anObject->GetName();
 135     if ( anObjName.isEmpty() )
 136       continue;
 137
 138     aNamesList.append( anObjName );
 139   }
 140
 141   QString aName;
 142
 143   if ( theIsTryToUsePurePrefix && !aNamesList.contains( thePrefix ) ) {
 144     aName = thePrefix;
 145   } else {
 146     int anId = 1;
 147     while( anId < aMaxNameId )
 148     {
 149       aName = QString( "%1_%2" ).arg( thePrefix ).arg( QString::number( anId++ ) );
 150
 151       // check that there are no other objects with the same name in the document
 152       if ( !aNamesList.contains( aName ) )
 153         break;
 154     }
 155   }
 156
 157   return aName;
 158 }
 159
 160 bool HYDROData_Tool::IsGeometryObject( const Handle(HYDROData_Entity)& theObject )
 161 {
 162   if ( theObject.IsNull() )
 163     return false;
 164
 165   return theObject->IsKind( STANDARD_TYPE(HYDROData_ArtificialObject) ) ||
 166          theObject->IsKind( STANDARD_TYPE(HYDROData_NaturalObject) );
 167 }
 168
 169 void HYDROData_Tool::UpdateChildObjectName( const QString&                  theOldStr,
 170                                             const QString&                  theNewStr,
 171                                             const Handle(HYDROData_Entity)& theObject )
 172 {
 173   if ( theObject.IsNull() )
 174     return;
 175
 176   QString anObjName = theObject->GetName();
 177   if ( theOldStr.isEmpty() )
 178   {
 179     while ( anObjName.startsWith( '_' ) )
 180       anObjName.remove( 0, 1 );
 181
 182     anObjName.prepend( theNewStr + "_" );
 183   }
 184   else if ( anObjName.startsWith( theOldStr ) )
 185   {
 186     anObjName.replace( 0, theOldStr.length(), theNewStr );
 187   }
 188   else
 189     return;
 190
 191   theObject->SetName( anObjName );
 192 }
 193
 194 QString HYDROData_Tool::GenerateNameForPython( const MapOfTreatedObjects& theTreatedObjects,
 195                                                const QString&             thePrefix )
 196 {
 197   QString aName = thePrefix;
 198   if ( !theTreatedObjects.contains( aName ) )
 199     return aName;
 200
 201   int anId = 1;
 202   while( anId < aMaxNameId )
 203   {
 204     aName = QString( "%1_%2" ).arg( thePrefix ).arg( QString::number( anId++ ) );
 205
 206     // check that there are no other objects with the same name
 207     if ( !theTreatedObjects.contains( aName ) )
 208       break;
 209   }
 210
 211   return aName;
 212 }
 213 //======================================================================================================
 214 TopAbs_State HYDROData_Tool::ComputePointState( const gp_XY& theXY, const TopoDS_Face& theFace )
 215 {
 216   TopAbs_State aState(TopAbs_UNKNOWN);
 217   if(theFace.IsNull()) return aState;
 218   Standard_Real aTol = BRep_Tool::Tolerance(theFace);
 219   BRepAdaptor_Surface Ads ( theFace, Standard_False );
 220   Standard_Real toluv = Min ( Ads.UResolution(aTol), Ads.VResolution(aTol) );
 221   if (toluv < 0.01)
 222     toluv = 0.01; // there is no need to be more precise than 1cm a any case !
 223                   // another solution could be to compute a tolerance related to the distance between the border nodes
 224   const gp_Pln& aPlane = Ads.Surface().Plane();
 225   gp_Pnt aPnt(theXY.X(), theXY.Y(), 0.);
 226   Standard_Real aU1, aV1;
 227   ElSLib::Parameters(aPlane,aPnt, aU1, aV1);
 228   BRepTopAdaptor_FClass2d aClassifier( theFace, toluv );
 229   aState = aClassifier.Perform( gp_Pnt2d(aU1, aV1), Standard_False );
 230   return aState;
 231 }
 232
 233 double HYDROData_Tool::GetAltitudeForEdge( const TopoDS_Edge& theEdge,
 234                                            const gp_XY& thePoint,
 235                                            double theParameterTolerance,
 236                                            double theSquareDistanceTolerance,
 237                                            double theInvalidAltitude )
 238 {
 239   double aFirst, aLast;
 240   Handle(Geom_Curve) aCurve = BRep_Tool::Curve( theEdge, aFirst, aLast );
 241   if( aCurve.IsNull() )
 242     return theInvalidAltitude;
 243
 244   gp_Pnt aFirstPnt, aLastPnt;
 245
 246   aCurve->D0( aFirst, aFirstPnt );
 247   aCurve->D0( aLast, aLastPnt );
 248
 249   gp_Pnt2d aFirstPnt2d( aFirstPnt.X(), aFirstPnt.Y() );
 250   gp_Pnt2d aLastPnt2d( aLastPnt.X(), aLastPnt.Y() );
 251
 252   double aFirstDist = 0;
 253   double aLastDist = aFirstPnt2d.SquareDistance( aLastPnt2d );
 254   double aNecDist = aFirstPnt2d.SquareDistance( thePoint );
 255
 256   while( fabs( aLast - aFirst ) > theParameterTolerance )
 257   {
 258     double aMid = ( aFirst + aLast ) / 2;
 259     gp_Pnt aMidPnt;
 260     aCurve->D0( aMid, aMidPnt );
 261     double aDist = aFirstPnt2d.SquareDistance( gp_Pnt2d( aMidPnt.X(), aMidPnt.Y() ) );
 262
 263     if( aDist < aNecDist )
 264       aFirst = aMid;
 265     else
 266       aLast = aMid;
 267   }
 268
 269   double aMid = ( aFirst + aLast ) / 2;
 270   gp_Pnt aMidPnt;
 271   aCurve->D0( aMid, aMidPnt );
 272
 273   gp_Pnt2d aMidPnt2d( aMidPnt.X(), aMidPnt.Y() );
 274   if( aMidPnt2d.SquareDistance( thePoint ) < theSquareDistanceTolerance )
 275     return aMidPnt.Z();
 276   else
 277     return theInvalidAltitude;
 278 }
 279
 280 double HYDROData_Tool::GetAltitudeForWire( const TopoDS_Wire& theWire,
 281                                            const gp_XY& thePoint,
 282                                            double theParameterTolerance,
 283                                            double theSquareDistanceTolerance,
 284                                            double theInvalidAltitude )
 285 {
 286   TopExp_Explorer anExp( theWire, TopAbs_EDGE );
 287   for( ; anExp.More(); anExp.Next() )
 288   {
 289     double anAltitude = GetAltitudeForEdge( TopoDS::Edge( anExp.Current() ), thePoint,
 290       theParameterTolerance, theSquareDistanceTolerance, theInvalidAltitude );
 291     if( anAltitude != theInvalidAltitude )
 292       return anAltitude;
 293   }
 294   return theInvalidAltitude;
 295 }
 296
 297 TopoDS_Shape HYDROData_Tool::getFirstShapeFromGroup( const HYDROData_SequenceOfObjects& theGroups,
 298                                                      const int                          theGroupId )
 299 {
 300   TopoDS_Shape aResShape;
 301   if ( theGroupId < 1 || theGroupId > theGroups.Length() )
 302     return aResShape;
 303
 304   Handle(HYDROData_ShapesGroup) aGroup =
 305     Handle(HYDROData_ShapesGroup)::DownCast( theGroups.Value( theGroupId ) );
 306   if ( aGroup.IsNull() )
 307     return aResShape;
 308
 309   TopTools_SequenceOfShape aGroupShapes;
 310   aGroup->GetShapes( aGroupShapes );
 311
 312   if ( !aGroupShapes.IsEmpty() )
 313     aResShape = aGroupShapes.First();
 314
 315   return aResShape;
 316 }
 317
 318 TCollection_ExtendedString HYDROData_Tool::toExtString( const QString& theStr )
 319 {
 320   TCollection_ExtendedString aRes;
 321   if( !theStr.isEmpty() )
 322   {
 323           Standard_ExtString extStr = new Standard_ExtCharacter[ ( theStr.length() + 1 ) * 2 ];
 324           memcpy( (void*)extStr, theStr.unicode(), theStr.length() * 2 );
 325           ((short*)extStr)[theStr.length()] = '\0';
 326     aRes = TCollection_ExtendedString( extStr );
 327           delete [] extStr;
 328   }
 329   return aRes;
 330 }
 331
 332 QString HYDROData_Tool::toQString( const TCollection_ExtendedString& theStr )
 333 {
 334   return QString( (QChar*)theStr.ToExtString(), theStr.Length() );
 335 }
 336
 337 Quantity_Color HYDROData_Tool::toOccColor( const QColor& theColor )
 338 {
 339   double r = theColor.red() / 255.0;
 340   double g = theColor.green() / 255.0;
 341   double b = theColor.blue() / 255.0;
 342
 343   return Quantity_Color( r, g, b, Quantity_TOC_RGB );
 344 }
 345
 346 QColor HYDROData_Tool::toQtColor( const Quantity_Color& theColor )
 347 {
 348   int r = 255 * theColor.Red();
 349   int g = 255 * theColor.Green();
 350   int b = 255 * theColor.Blue();
 351   return QColor( r, g, b );
 352 }
 353
 354 QColor HYDROData_Tool::GenerateRandColor()
 355 {
 356   float aHue = ( rand()%1000 ) * 0.001f;
 357
 358   QColor aColor;
 359   aColor.setHsl( (int)(aHue*255.), 200, 128 );
 360   int r = aColor.red();
 361   int g = aColor.green();
 362   int b = aColor.blue();
 363   return ( aColor.isValid() ? aColor : Qt::darkBlue );
 364 }
 365
 366 void HYDROData_Tool::GenerateRepeatableRandColors(int nbColorsToGen, QVector<QColor>& theColors)
 367 {
 368   for (int i = 1; i <= nbColorsToGen; i++)
 369     theColors.append(HYDROData_Tool::GenerateRandColor());
 370 }
 371
 372 bool HYDROData_Tool::GenerateNonRepeatableRandColors(int nbColorsToGen, QVector<QColor>& theColors)
 373 {
 374   if (nbColorsToGen > aMaxColorNb)
 375     return false;
 376   QSet<int> Codes;
 377   float aHue;
 378   theColors.clear();
 379   do
 380   {
 381     QColor aColor;
 382     int H = rand()%360;
 383     int S = rand()%256;
 384     int Code = H*256+S;
 385     if (Codes.contains(Code))
 386       continue;
 387     aColor.setHsl( H, S, 128 );
 388     if (aColor.isValid())
 389     {
 390       theColors.append(aColor);
 391       Codes.insert(Code);
 392     }
 393   } while (theColors.size() <= nbColorsToGen);
 394   return true;
 395 }
 396
 397
 398 bool HYDROData_Tool::IsNan( double theValue )
 399 {
 400 #ifdef WIN32
 401   return _isnan( theValue );
 402 #else
 403   return isnan( theValue );
 404 #endif
 405 }
 406
 407 bool HYDROData_Tool::IsInf( double theValue )
 408 {
 409 #ifdef WIN32
 410   return (!_finite( theValue  ) );
 411 #else
 412   return isinf( theValue );
 413 #endif
 414 }
 415
 416 static void MakeShellG(const NCollection_Map<TopoDS_Face, TopTools_ShapeMapHasher>& FG,
 417   TopoDS_Shape& outSh)
 418 {
 419   BRep_Builder bb;
 420   NCollection_Map<TopoDS_Face, TopTools_ShapeMapHasher>::Iterator itFG(FG);
 421   if (FG.Extent() > 1)
 422   {
 423     //face nb > 1 => make shell
 424     TopoDS_Shell outShell;
 425     bb.MakeShell(outShell);
 426     for (;itFG.More();itFG.Next())
 427       bb.Add(outShell, itFG.Value());
 428     outSh = outShell;
 429   }
 430   else if (FG.Extent() == 1)
 431   {
 432     outSh = itFG.Value(); //one face
 433   }
 434 }
 435
 436 TopoDS_Shape HYDROData_Tool::RebuildCmp(const TopoDS_Shape& in)
 437 {
 438   TopTools_IndexedDataMapOfShapeListOfShape mE2LF;
 439   TopExp::MapShapesAndAncestors(in, TopAbs_EDGE, TopAbs_FACE, mE2LF);
 440   if (mE2LF.IsEmpty())
 441     return TopoDS_Shape();
 442   NCollection_Map<TopoDS_Face, TopTools_ShapeMapHasher> dfm;
 443   //TopExp::MapShapes(aFuseShape, TopAbs_FACE, dfm);
 444   TopExp_Explorer expf(in, TopAbs_FACE);
 445   for (;expf.More(); expf.Next())
 446     dfm.Add(TopoDS::Face(expf.Current()));
 447
 448   int nbF = dfm.Extent();
 449   TopExp_Explorer exp_f(in, TopAbs_FACE);
 450   const TopoDS_Face& FF = TopoDS::Face(exp_f.Current());
 451   NCollection_List<TopoDS_Face> CurrFS;
 452   NCollection_List<TopoDS_Face> NeighFS;
 453   NCollection_Map<TopoDS_Face, TopTools_ShapeMapHasher> PrF;
 454   CurrFS.Append(FF);
 455   NCollection_List<NCollection_Map<TopoDS_Face, TopTools_ShapeMapHasher>> GL_F;
 456   NCollection_Map<TopoDS_Face, TopTools_ShapeMapHasher> OneGr;
 457   bool end = false;
 458   while (!end)
 459   {
 460     NCollection_List<TopoDS_Face>::Iterator it_currfs(CurrFS);
 461     NeighFS.Clear();
 462     for (;it_currfs.More();it_currfs.Next())
 463     {
 464       const TopoDS_Face& CF = it_currfs.Value();
 465       TopExp_Explorer exp_edge(CF, TopAbs_EDGE);
 466       for (;exp_edge.More();exp_edge.Next())
 467       {
 468         const TopoDS_Shape& CE = exp_edge.Current();
 469         const TopTools_ListOfShape& lsf = mE2LF.FindFromKey(CE);
 470         TopTools_ListIteratorOfListOfShape ls_it(lsf); //always one face (since all faces are planar)
 471         for (;ls_it.More();ls_it.Next())
 472         {
 473           const TopoDS_Face& F = TopoDS::Face(ls_it.Value());
 474           if (F.IsSame(CF))
 475             continue;
 476           if (!PrF.Contains(F))
 477           {
 478             OneGr.Add(F);
 479             NeighFS.Append(F);
 480             PrF.Add(F);
 481           }
 482         }
 483       }
 484       OneGr.Add(CF);
 485       PrF.Add(CF);
 486     }
 487     if (NeighFS.IsEmpty())
 488     {
 489       GL_F.Append(OneGr);
 490       OneGr.Clear();
 491       dfm.Subtract(PrF);
 492       if (dfm.IsEmpty())
 493         end = true;
 494       else
 495       {
 496         NCollection_Map<TopoDS_Face, TopTools_ShapeMapHasher>::Iterator itDm(dfm);
 497         const TopoDS_Face& nsh = itDm.Key();
 498         NeighFS.Append(nsh);
 499       }
 500     }
 501     CurrFS = NeighFS;
 502   }
 503
 504   TopoDS_Shape sh;
 505
 506   if (GL_F.Extent() > 1)
 507   {
 508     TopoDS_Compound cmp;
 509     NCollection_List<NCollection_Map<TopoDS_Face, TopTools_ShapeMapHasher>>::Iterator itGL_F(GL_F);
 510     BRep_Builder bb;
 511     bb.MakeCompound(cmp);
 512     for (;itGL_F.More();itGL_F.Next())
 513     {
 514       MakeShellG(itGL_F.Value(), sh);
 515       if (!sh.IsNull())
 516         bb.Add(cmp, sh);
 517     }
 518     return  cmp;
 519   }
 520   else if (GL_F.Extent() == 1)
 521   {
 522     MakeShellG(GL_F.First(), sh);
 523     return sh;
 524   }
 525
 526 }
 527
 528 TopoDS_Shape HYDROData_Tool::PolyXY2Face( const Handle(HYDROData_PolylineXY)& aPolyline )
 529 {
 530   //DEBTRACE("generateTopShape");
 531   TopoDS_Face aResultFace = TopoDS_Face(); // --- result: default = null face
 532
 533   if (!aPolyline.IsNull())
 534     {
 535       TopoDS_Shape aPolylineShape = aPolyline->GetShape();
 536 #ifdef DEB_IMZ
 537       std::string brepName = "imz.brep";
 538       BRepTools::Write(aPolylineShape, brepName.c_str());
 539 #endif
 540       TopTools_ListOfShape aWiresList;
 541
 542       if (!aPolylineShape.IsNull() && aPolylineShape.ShapeType() == TopAbs_WIRE)
 543         {
 544           // --- only one wire: try to make a face
 545           //DEBTRACE("one wire: try to build a face");
 546           const TopoDS_Wire& aPolylineWire = TopoDS::Wire(aPolylineShape);
 547           if (!aPolylineWire.IsNull())
 548             {
 549               BRepBuilderAPI_MakeFace aMakeFace(aPolylineWire, Standard_True);
 550               aMakeFace.Build();
 551               if (aMakeFace.IsDone())
 552                 {
 553                   //DEBTRACE(" a face with the only wire given");
 554                   aResultFace = aMakeFace.Face();
 555                 }
 556             }
 557         }
 558       else
 559         {
 560           // --- a list of wires ? inventory of wires and edges
 561           Handle(TopTools_HSequenceOfShape) aSeqWires = new TopTools_HSequenceOfShape;
 562           Handle(TopTools_HSequenceOfShape) aSeqEdges = new TopTools_HSequenceOfShape;
 563           TopExp_Explorer anExp(aPolylineShape, TopAbs_WIRE);
 564           //DEBTRACE("list of wires ?");
 565           for (; anExp.More(); anExp.Next())
 566             {
 567               if (!anExp.Current().IsNull())
 568                 {
 569                   const TopoDS_Wire& aWire = TopoDS::Wire(anExp.Current());
 570                   aWiresList.Append(aWire);
 571                   //DEBTRACE("  append wire");
 572                   TopExp_Explorer it2(aWire, TopAbs_EDGE);
 573                   for (; it2.More(); it2.Next())
 574                     aSeqEdges->Append(it2.Current());
 575                 }
 576             }
 577           if (aWiresList.IsEmpty())
 578             return aResultFace; // --- no wires: null result
 579
 580           if (aSeqEdges->Length() > 1)
 581             {
 582               //DEBTRACE("try to connect all the edges together, build a unique wire and a face");
 583               // --- try to create one wire by connecting edges with a distance tolerance (no necessity of sharing points)
 584               ShapeAnalysis_FreeBounds::ConnectEdgesToWires(aSeqEdges, 1E-5, Standard_False, aSeqWires);
 585
 586               if (aSeqWires->Length() == 1)
 587                 {
 588                   // --- one wire: try to make a face
 589                   const TopoDS_Wire& aPolylineWire = TopoDS::Wire(aSeqWires->Value(1));
 590                   if (!aPolylineWire.IsNull())
 591                     {
 592                       BRepBuilderAPI_MakeFace aMakeFace(aPolylineWire, Standard_True);
 593                       aMakeFace.Build();
 594                       if (aMakeFace.IsDone())
 595                         {
 596                           //DEBTRACE("  a face from all the wires connected");
 597                           aResultFace = aMakeFace.Face();
 598                         }
 599                     }
 600                 }
 601             }
 602
 603           if (aResultFace.IsNull())
 604             {
 605               //DEBTRACE("try to make a face with the first wire of the list and other wires as restrictions");
 606               // --- try to make a face with the first wire of the list and other wires as restrictions
 607               BRepAlgo_FaceRestrictor aFR;
 608               TopoDS_Face aRefFace;
 609               TopoDS_Shape aS = aWiresList.First();
 610               BRepBuilderAPI_MakeFace aMakeFace(TopoDS::Wire(aWiresList.First()), Standard_True);
 611               aMakeFace.Build();
 612               if (aMakeFace.IsDone())
 613                 {
 614                   //DEBTRACE("  a face with first wire");
 615                   aRefFace = aMakeFace.Face();
 616                 }
 617               if (!aRefFace.IsNull())
 618                 {
 619                   aFR.Init(aRefFace, Standard_False, Standard_True);
 620                   TopTools_ListIteratorOfListOfShape anIt(aWiresList);
 621                   for (; anIt.More(); anIt.Next())
 622                     {
 623                       TopoDS_Wire& aWire = TopoDS::Wire(anIt.Value());
 624                       if (aWire.IsNull())
 625                         continue;
 626                       aFR.Add(aWire);
 627                     }
 628                   aFR.Perform();
 629                   if (aFR.IsDone())
 630                     {
 631                       for (; aFR.More(); aFR.Next())
 632                         {
 633                           //DEBTRACE("  a restricted face");
 634                           aResultFace = aFR.Current();
 635                           break;
 636                         }
 637                     }
 638                 }
 639             }
 640         }
 641     }
 642
 643   if (aResultFace.IsNull())
 644     return aResultFace;
 645
 646   //DEBTRACE("check the face");
 647   BRepCheck_Analyzer anAnalyzer(aResultFace);
 648   if (anAnalyzer.IsValid() && aResultFace.ShapeType() == TopAbs_FACE)
 649   {
 650     //DEBTRACE("face OK");
 651     return aResultFace;
 652   }
 653   else
 654   {
 655     //DEBTRACE("bad face");
 656   }
 657   return TopoDS_Face();
 658 }
 659
 660 void HYDROData_Tool::SetSIProgress(const Handle(Message_ProgressIndicator)& thePI)
 661 {
 662   StricklerInterpolationProgress() = thePI;
 663 }
 664
 665 const Handle(Message_ProgressIndicator)& HYDROData_Tool::GetSIProgress()
 666 {
 667   return StricklerInterpolationProgress();
 668 }
 669
 670 Handle(Message_ProgressIndicator)& HYDROData_Tool::StricklerInterpolationProgress()
 671 {
 672   static Handle(Message_ProgressIndicator) aPI = NULL;
 673   return aPI;
 674 }
 675
 676 void HYDROData_Tool::SetZIProgress(const Handle(Message_ProgressIndicator)& thePI)
 677 {
 678   BathymetryInterpolationProgress() = thePI;
 679 }
 680
 681 const Handle(Message_ProgressIndicator)& HYDROData_Tool::GetZIProgress()
 682 {
 683   return BathymetryInterpolationProgress();
 684 }
 685
 686 Handle(Message_ProgressIndicator)& HYDROData_Tool::BathymetryInterpolationProgress()
 687 {
 688   static Handle(Message_ProgressIndicator) aPI = NULL;
 689   return aPI;
 690 }
 691
 692 void HYDROData_Tool::SetTriangulationStatus(const ExecStatus& theStatus)
 693 {
 694   myTriangulationStatus = theStatus;
 695 }
 696
 697 const HYDROData_Tool::ExecStatus& HYDROData_Tool::GetTriangulationStatus()
 698 {
 699   return myTriangulationStatus;
 700 }
 701
 702 static bool AddXYZ(bool bImportXY,
 703                    double x,
 704                    double y,
 705                    double z,
 706                    std::vector<gp_XYZ>& thePointsXYZ,
 707                    std::vector<gp_XY>& thePointsXY)
 708 {
 709   if (!bImportXY)
 710   {
 711     if ( HYDROData_Tool::IsNan( x ) || HYDROData_Tool::IsInf( x ) ||
 712       HYDROData_Tool::IsNan( y ) || HYDROData_Tool::IsInf( y ) ||
 713       HYDROData_Tool::IsNan( z ) || HYDROData_Tool::IsInf( z ) )
 714       return false;
 715
 716     if( thePointsXYZ.size()>=thePointsXYZ.capacity() )
 717       thePointsXYZ.reserve( thePointsXYZ.size()+BLOCK_SIZE );
 718
 719     thePointsXYZ.push_back(gp_XYZ(x,y,z));
 720   }
 721   else
 722   {
 723     if ( HYDROData_Tool::IsNan( x ) || HYDROData_Tool::IsInf( x ) ||
 724       HYDROData_Tool::IsNan( y ) || HYDROData_Tool::IsInf( y ) )
 725       return false;
 726
 727     if( thePointsXY.size()>=thePointsXY.capacity() )
 728       thePointsXY.reserve( thePointsXY.size()+BLOCK_SIZE );
 729
 730     thePointsXY.push_back(gp_XY(x,y));
 731   }
 732   return true;
 733 }
 734
 735 bool HYDROData_Tool::importFromXYZ( QString& theFileName,
 736                                     bool bImportXY,
 737                                     std::vector<gp_XYZ>& thePointsXYZ,
 738                                     std::vector<gp_XY>& thePointsXY)
 739 {
 740   QFile aFile( theFileName );
 741   if ( !aFile.exists() || !aFile.open( QIODevice::ReadOnly ) )
 742     return false;
 743
 744   QString aFileSuf = QFileInfo( aFile ).suffix().toLower();
 745
 746   double x,y;
 747   if ( aFileSuf == "xyz" )
 748   {
 749     double z;
 750     while ( !aFile.atEnd() )
 751     {
 752       std::string aLine = aFile.readLine().simplified().toStdString();
 753       if ( aLine.empty() )
 754         continue;
 755
 756       x = 0;
 757       y = 0;
 758       z = 0;
 759
 760       if( sscanf( aLine.c_str(), "%lf %lf %lf", &x, &y, &z )!=3 )
 761       {
 762         aFile.close();
 763         return false;
 764       }
 765
 766       if (!AddXYZ(bImportXY, x, y, z, thePointsXYZ, thePointsXY ))
 767       {
 768         aFile.close();
 769         return false;
 770       }
 771     }
 772   }
 773   else if (aFileSuf == "xy" )
 774   {
 775     while ( !aFile.atEnd() )
 776     {
 777       std::string aLine = aFile.readLine().simplified().toStdString();
 778       if ( aLine.empty() )
 779         continue;
 780
 781       x = 0;
 782       y = 0;
 783
 784       if( sscanf( aLine.c_str(), "%lf %lf", &x, &y )!=2 )
 785       {
 786         aFile.close();
 787         return false;
 788       }
 789
 790       if (!AddXYZ(true, x, y, 0, thePointsXYZ, thePointsXY ))
 791       {
 792         aFile.close();
 793         return false;
 794       }
 795     }
 796   }
 797
 798   aFile.close();
 799
 800   return true;
 801 }
 802
 803 bool HYDROData_Tool::importPolylineFromXYZ(QString aFileName, Handle(HYDROData_Document) theDocument,
 804   bool importXY, NCollection_Sequence<Handle(HYDROData_Entity)>& importedEntities)
 805 {
 806   if (importXY)
 807   {
 808     std::vector<gp_XY> aPoints2d;
 809     std::vector<gp_XYZ> aDPoints3d;
 810
 811     if (HYDROData_Tool::importFromXYZ(aFileName, importXY, aDPoints3d, aPoints2d))
 812     {
 813       QString basename = QFileInfo( aFileName ).baseName();
 814
 815       Handle(HYDROData_PolylineXY) aPolylineXY = Handle(HYDROData_PolylineXY)::DownCast( theDocument->CreateObject( KIND_POLYLINEXY ) );
 816       HYDROData_PolylineXY::SectionType aSectType = HYDROData_PolylineXY::SECTION_POLYLINE;
 817       bool IsClosed = false;
 818       if ((aPoints2d.front()-aPoints2d.back()).Modulus()<Precision::Confusion())
 819       {
 820         IsClosed = true;
 821         aPolylineXY->AddSection( TCollection_AsciiString("poly_section"), aSectType, true);
 822       }
 823       else
 824         aPolylineXY->AddSection( TCollection_AsciiString("poly_section"), aSectType, false);
 825
 826       int n = aPoints2d.size();
 827       if (IsClosed)
 828         n--;
 829
 830       for ( int i = 0; i < n; i++ )
 831       {
 832         gp_XY aSectPoint = aPoints2d[i];
 833         theDocument->Transform(aSectPoint, true);
 834         aPolylineXY->AddPoint( 0, aSectPoint );
 835       }
 836
 837       aPolylineXY->SetWireColor( HYDROData_PolylineXY::DefaultWireColor() );
 838       aPolylineXY->SetName( basename + "_PolyXY_" );
 839       aPolylineXY->Update();
 840       importedEntities.Append(aPolylineXY);
 841       return true;
 842     }
 843     else
 844       return false;
 845   }
 846   else //xyz
 847   {
 848     std::vector<gp_XY> aDPoints2d;
 849     std::vector<gp_XYZ> aPoints3d;
 850     if (HYDROData_Tool::importFromXYZ(aFileName, false, aPoints3d, aDPoints2d))
 851     {
 852       QString basename = QFileInfo( aFileName ).baseName();
 853       Handle(HYDROData_PolylineXY) aPolylineXY = Handle(HYDROData_PolylineXY)::DownCast( theDocument->CreateObject( KIND_POLYLINEXY ) );
 854       Handle(HYDROData_Polyline3D) aPolylineObj = Handle(HYDROData_Polyline3D)::DownCast( theDocument->CreateObject( KIND_POLYLINE ) );
 855       Handle(HYDROData_Bathymetry) aBath = Handle(HYDROData_Bathymetry)::DownCast( theDocument->CreateObject( KIND_BATHYMETRY ) );
 856       HYDROData_Bathymetry::AltitudePoints aAPoints;
 857       HYDROData_PolylineXY::SectionType aSectType = HYDROData_PolylineXY::SECTION_POLYLINE;
 858       bool IsClosed = false;
 859       if ((aPoints3d.front()-aPoints3d.back()).Modulus()<Precision::Confusion())
 860       {
 861         IsClosed = true;
 862         aPolylineXY->AddSection( TCollection_AsciiString("poly_section"), aSectType, true);
 863       }
 864       else
 865         aPolylineXY->AddSection( TCollection_AsciiString("poly_section"), aSectType, false);
 866
 867       int n = aPoints3d.size();
 868       if (IsClosed)
 869         n--;
 870
 871       for ( int i = 0; i < n; i++ )
 872       {
 873         gp_XY aSectPoint(aPoints3d[i].X(), aPoints3d[i].Y());
 874         theDocument->Transform(aSectPoint, true);
 875         aPolylineXY->AddPoint( 0, aSectPoint );
 876         HYDROData_Bathymetry::AltitudePoint p;
 877         p.X = aSectPoint.X();
 878         p.Y = aSectPoint.Y();
 879         p.Z = aPoints3d[i].Z();
 880         aAPoints.push_back(p);
 881       }
 882
 883       QString aBathName = basename + "_bath_";
 884       QString aPolyXYName = basename + "_polyXY_";
 885       QString aPoly3DName = basename + "_poly3D_";
 886
 887       aPolylineXY->SetName( aPolyXYName );
 888       aPolylineXY->SetWireColor(HYDROData_PolylineXY::DefaultWireColor());
 889       aPolylineXY->Update();
 890
 891       aBath->SetAltitudePoints(aAPoints);
 892       aBath->SetName( aBathName );
 893
 894       aPolylineObj->SetPolylineXY (aPolylineXY, false);
 895       aPolylineObj->SetAltitudeObject(aBath);
 896
 897       aPolylineObj->SetBorderColor( aPolylineObj->DefaultBorderColor() );
 898       aPolylineObj->SetName( aPoly3DName );
 899
 900       aPolylineObj->Update();
 901       importedEntities.Append(aPolylineXY);
 902       importedEntities.Append(aPolylineObj);
 903       return true;
 904     }
 905     else
 906       return false;
 907   }
 908 }
 909
 910
 911
 912
 913 Handle(Geom2d_Curve) HYDROData_Tool::BRepAdaptorTo2DCurve( const BRepAdaptor_Curve& ad )
 914 {
 915   if( ad.GetType() == GeomAbs_Line)
 916   {
 917     double f = ad.FirstParameter();
 918     double l = ad.LastParameter();
 919     return new Geom2d_TrimmedCurve( GeomAPI::To2d(ad.Curve().Curve(), Geom_Plane(gp::XOY()).Pln()), f, l );
 920   }
 921
 922   if( ad.GetType() == GeomAbs_BSplineCurve )
 923   {
 924     //just ignore Z coord, without projecting
 925     Handle(Geom_BSplineCurve) aSpline = ad.Curve().BSpline();
 926     if (aSpline.IsNull())
 927       return Handle(Geom2d_Curve)();
 928     TColgp_Array1OfPnt poles3d = aSpline->Poles();
 929     TColgp_HArray1OfPnt2d poles2d(poles3d.Lower(), poles3d.Upper());
 930     for (int i=poles3d.Lower(); i<=poles3d.Upper();i++)
 931     {
 932       gp_XY p2d(poles3d(i).X(), poles3d(i).Y());
 933       poles2d(i).SetXY(p2d);
 934     }
 935     const TColStd_Array1OfReal& knots = aSpline->Knots();
 936     const TColStd_Array1OfInteger& multiplicities = aSpline->Multiplicities();
 937     int aDegree = aSpline->Degree();
 938     return new Geom2d_BSplineCurve( poles2d, knots, multiplicities, aDegree );
 939   }
 940   return Handle(Geom2d_Curve)();
 941 }
 942
 943 std::ostream& operator<<( std::ostream& theStream, const QString& theText )
 944 {
 945   theStream << theText.toStdString();
 946   return theStream;
 947 }
 948
 949 std::ostream& operator<<( std::ostream& theStream, const QColor& theColor )
 950 {
 951   theStream << "[" << theColor.red() << ", " << theColor.green() << ", " << theColor.blue() << "]";
 952   return theStream;
 953 }
 954
 955 std::ostream& operator<<( std::ostream& theStream, const TopoDS_Shape& theShape )
 956 {
 957   theStream << "[" << theShape.TShape().operator->() << "]";
 958   return theStream;
 959 }
 960
 961 std::ostream& operator<<( std::ostream& theStream, const TopoDS_Face& theFace )
 962 {
 963   theStream << "[" << theFace.TShape().operator->() << "]";
 964   return theStream;
 965 }
 966
 967 std::ostream& operator<<( std::ostream& theStream, const gp_XY& theXY )
 968 {
 969   theStream << "(" << theXY.X() << "; " << theXY.Y() << ")";
 970   return theStream;
 971 }
 972
 973 bool operator == ( const gp_XY& thePoint1, const gp_XY& thePoint2 )
 974 {
 975   const double EPS = 1E-3;
 976   return
 977     fabs( thePoint1.X() - thePoint2.X() ) < EPS &&
 978     fabs( thePoint1.Y() - thePoint2.Y() ) < EPS;
 979
 980 }