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.
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.
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
16 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 #include <HYDROData_LandCoverMap.h>
20 #include <HYDROData_Object.h>
21 #include <HYDROData_PolylineXY.h>
22 #include <HYDROData_Tool.h>
24 #include <BOPAlgo_BOP.hxx>
25 #include <BOPAlgo_Builder.hxx>
26 #include <BOPAlgo_PaveFiller.hxx>
27 #include <BOPCol_ListOfShape.hxx>
28 #include <BRep_Builder.hxx>
29 #include <BRepAlgoAPI_Fuse.hxx>
30 #include <BRepBuilderAPI_MakeFace.hxx>
31 #include <NCollection_IndexedMap.hxx>
33 #include <TopoDS_Compound.hxx>
34 #include <TopoDS_Edge.hxx>
35 #include <TopoDS_Face.hxx>
36 #include <TopoDS_Iterator.hxx>
37 #include <TopoDS_Shell.hxx>
38 #include <TopExp_Explorer.hxx>
39 #include <TopTools_ListIteratorOfListOfShape.hxx>
43 #include <QTextStream>
45 const char TELEMAC_FORMAT = 'f';
46 const int TELEMAC_PRECISION = 3;
49 IMPLEMENT_STANDARD_HANDLE(HYDROData_LandCoverMap, HYDROData_Entity)
50 IMPLEMENT_STANDARD_RTTIEXT(HYDROData_LandCoverMap, HYDROData_Entity)
52 class HYDROData_MapOfFaceToStricklerType : public NCollection_IndexedDataMap<TopoDS_Face, QString>
58 @param theMap the land cover map to iterate through
60 HYDROData_LandCoverMap::Iterator::Iterator( const HYDROData_LandCoverMap& theMap )
65 HYDROData_LandCoverMap::Iterator::Iterator( const Handle( HYDROData_LandCoverMap )& theMap )
77 Initialize the iterator
78 @param theMap the land cover map to iterate through
80 void HYDROData_LandCoverMap::Iterator::Init( const HYDROData_LandCoverMap& theMap )
82 TopoDS_Shape aShape = theMap.GetShape();
86 myIterator = new TopoDS_Iterator( aShape );
88 theMap.myLab.FindChild( DataTag_Types ).FindAttribute( TDataStd_ExtStringArray::GetID(), myArray );
94 HYDROData_LandCoverMap::Iterator::~Iterator()
100 Return if the iterator has more elements
101 @return if the iterator has more elements
103 bool HYDROData_LandCoverMap::Iterator::More() const
105 return !myArray.IsNull() && myIterator && myIterator->More();
109 Move iterator to the next element
111 void HYDROData_LandCoverMap::Iterator::Next()
121 Get the current land cover (face)
122 @return the land cover's face
124 TopoDS_Face HYDROData_LandCoverMap::Iterator::Face() const
127 return TopoDS::Face( myIterator->Value() );
129 return TopoDS_Face();
133 Get the current land cover's Strickler type
134 @return the land cover's Strickler type
136 QString HYDROData_LandCoverMap::Iterator::StricklerType() const
138 if( myArray.IsNull() || myIndex < myArray->Lower() || myIndex > myArray->Upper() )
141 return HYDROData_Tool::toQString( myArray->Value( myIndex ) );
147 HYDROData_LandCoverMap::HYDROData_LandCoverMap()
148 : HYDROData_Entity( Geom_No )
155 HYDROData_LandCoverMap::~HYDROData_LandCoverMap()
161 @return object's kind
163 const ObjectKind HYDROData_LandCoverMap::GetKind() const
165 return KIND_LAND_COVER_MAP;
169 Import the land cover map from QGIS
170 @param theFileName the name of file
171 @return if the import is successful
173 bool HYDROData_LandCoverMap::ImportQGIS( const QString& theFileName )
180 Export the land cover map to QGIS
181 @param theFileName the name of file
182 @return if the export is successful
184 bool HYDROData_LandCoverMap::ExportQGIS( const QString& theFileName ) const
190 void EdgeDiscretization( const TopoDS_Edge& theEdge,
191 NCollection_IndexedMap<gp_Pnt>& theVerticesMap,
192 QList<int>& theVerticesIds )
198 Export the land cover map for the solver (Telemac)
199 @param theFileName the name of file
200 @return if the export is successful
202 bool HYDROData_LandCoverMap::ExportTelemac( const QString& theFileName ) const
204 TopoDS_Shell aShell; //TODO: unite all the faces of land covers into the shell
206 NCollection_IndexedMap<gp_Pnt> aVerticesMap;
207 NCollection_IndexedDataMap< TopoDS_Edge, QList<int> > anEdgesMap;
208 NCollection_IndexedDataMap< TopoDS_Face, QList<int> > aFacesMap;
210 // add into the map all edges existing in the shell
211 TopExp_Explorer anExp1( aShell, TopAbs_EDGE );
212 for( ; anExp1.More(); anExp1.Next() )
214 TopoDS_Edge anEdge = TopoDS::Edge( anExp1.Current() );
215 QList<int> aVerticesIdsList;
216 EdgeDiscretization( anEdge, aVerticesMap, aVerticesIdsList );
217 anEdgesMap.Add( anEdge, aVerticesIdsList );
220 // add into the map all faces existing in the shell and correspondence between face and edges ids
221 TopExp_Explorer anExp2( aShell, TopAbs_FACE );
222 for( ; anExp2.More(); anExp2.Next() )
224 TopoDS_Face aFace = TopoDS::Face( anExp2.Current() );
225 TopExp_Explorer anExp3( aFace, TopAbs_EDGE );
226 QList<int> anEdgesIdsList;
227 for( ; anExp3.More(); anExp3.Next() )
229 TopoDS_Edge anEdge = TopoDS::Edge( anExp3.Current() );
230 int anEdgeId = anEdgesMap.FindIndex( anEdge );
231 anEdgesIdsList.append( anEdgeId );
233 aFacesMap.Add( aFace, anEdgesIdsList );
236 QFile aFile( theFileName );
237 if( !aFile.open( QFile::WriteOnly | QFile::Text ) )
240 QTextStream aStream( &aFile );
241 aStream << "# nodes\n";
242 NCollection_IndexedMap<gp_Pnt>::Iterator anIt1( aVerticesMap );
243 for( ; anIt1.More(); anIt1.Next() )
245 gp_Pnt aPnt = anIt1.Value();
246 aStream << QString::number( aPnt.X(), TELEMAC_FORMAT, TELEMAC_PRECISION );
248 aStream << QString::number( aPnt.Y(), TELEMAC_FORMAT, TELEMAC_PRECISION );
250 aStream << QString::number( aPnt.Z(), TELEMAC_FORMAT, TELEMAC_PRECISION );
255 aStream << "# edges\n";
256 NCollection_IndexedDataMap< TopoDS_Edge, QList<int> >::Iterator anIt2( anEdgesMap );
257 for( ; anIt2.More(); anIt2.Next() )
259 QList<int> aVerticesIds = anIt2.Value();
260 foreach( int anId, aVerticesIds )
261 aStream << anId << " ";
266 aStream << "# faces\n";
267 NCollection_IndexedDataMap< TopoDS_Face, QList<int> >::Iterator anIt3( aFacesMap );
268 for( ; anIt3.More(); anIt3.Next() )
270 QList<int> anEdgesIds = anIt3.Value();
271 foreach( int anId, anEdgesIds )
272 aStream << anId << " ";
282 Add a new object as land cover
283 @param theObject the object to add as land cover
284 @param theType the Strickler type for the new land cover
285 @return if the addition is successful
287 bool HYDROData_LandCoverMap::Add( const Handle( HYDROData_Object )& theObject, const QString& theType )
289 if( theObject.IsNull() )
292 TopoDS_Shape aShape = theObject->GetTopShape();
293 if( aShape.ShapeType()!=TopAbs_FACE )
296 TopoDS_Face aFace = TopoDS::Face( aShape );
297 return LocalPartition( aFace, theType );
301 Add a new polyline as land cover
302 @param thePolyline the polyline to add as land cover
303 @param theType the Strickler type for the new land cover
304 @return if the addition is successful
306 bool HYDROData_LandCoverMap::Add( const Handle( HYDROData_PolylineXY )& thePolyline, const QString& theType )
308 if( thePolyline.IsNull() )
311 TopoDS_Shape aShape = thePolyline->GetShape();
312 if( aShape.ShapeType()!=TopAbs_WIRE )
315 TopoDS_Wire aWire = TopoDS::Wire( aShape );
316 if( !aWire.Closed() )
319 TopoDS_Face aFace = BRepBuilderAPI_MakeFace( aWire, Standard_True ).Face();
320 return LocalPartition( aFace, theType );
324 Remove the given face from land cover map
325 @param theFace the face to be removed
326 @return if the removing is successful
328 bool HYDROData_LandCoverMap::Remove( const TopoDS_Face& theFace )
330 TopTools_ListOfShape aList;
331 aList.Append( theFace );
332 return Remove( aList );
336 Remove the given faces from land cover map
337 @param theFacesToRemove the face list to be removed
338 @return if the removing is successful
340 bool HYDROData_LandCoverMap::Remove( const TopTools_ListOfShape& theFacesToRemove )
342 HYDROData_MapOfFaceToStricklerType aFacesToRemove, aNewFaces;
343 TopTools_ListIteratorOfListOfShape aFIt( theFacesToRemove );
344 for( ; aFIt.More(); aFIt.Next() )
346 TopoDS_Shape aShape = aFIt.Value();
347 if( aShape.ShapeType()==TopAbs_FACE )
348 aFacesToRemove.Add( TopoDS::Face( aShape ), "" );
351 Iterator anIt( *this );
352 for( ; anIt.More(); anIt.Next() )
353 if( !aFacesToRemove.Contains( anIt.Face() ) )
354 aNewFaces.Add( anIt.Face(), anIt.StricklerType() );
356 StoreLandCovers( aNewFaces );
361 Split the land cover map by the given polyline
362 @param thePolyline the tool polyline to split the land cover map
363 @return if the removing is successful
365 bool HYDROData_LandCoverMap::Split( const Handle( HYDROData_PolylineXY )& thePolyline )
367 if( thePolyline.IsNull() )
370 TopoDS_Shape aShape = thePolyline->GetShape();
371 return LocalPartition( aShape, "" );
375 Merge the given faces in the land cover
376 @param theFaces the faces to merge in the land cover map
377 @param theType the Strickler type for the merged land cover
378 @return if the merge is successful
380 bool HYDROData_LandCoverMap::Merge( const TopTools_ListOfShape& theFaces, const QString& theType )
382 // 1. to remove the merged faces from the current map
385 // 2. to fuse the faces into the new face
386 BOPAlgo_PaveFiller aPF;
387 aPF.SetArguments( theFaces );
388 aPF.SetFuzzyValue( 1E-2 );
389 aPF.SetRunParallel( Standard_False );
391 int iErr = aPF.ErrorStatus();
396 aBOP.SetArguments( theFaces );
397 aBOP.SetOperation( BOPAlgo_FUSE );
398 aBOP.SetRunParallel( Standard_False );
399 aBOP.PerformWithFiller(aPF);
400 iErr = aBOP.ErrorStatus();
404 TopoDS_Shape aMergedShape = aBOP.Shape();
405 if( aMergedShape.ShapeType()!=TopAbs_FACE )
408 // 3. to add the face into the map
409 return LocalPartition( TopoDS::Face( aMergedShape ), theType );
413 Get the shape of the land cover map
415 TopoDS_Shape HYDROData_LandCoverMap::GetShape() const
417 return HYDROData_Entity::GetShape( DataTag_Shape );
421 Set the shape of the land cover map
422 @param theShape the new shape for the land cover map
424 void HYDROData_LandCoverMap::SetShape( const TopoDS_Shape& theShape )
426 HYDROData_Entity::SetShape( DataTag_Shape, theShape );
430 Perform the local partition algorithm on the land cover
431 @param theNewShape the new shape to add into the land cover
432 @param theNewType the new Strickler type for the new land cover
433 @return if the local partition is successful
435 bool HYDROData_LandCoverMap::LocalPartition( const TopoDS_Shape& theNewShape, const QString& theNewType )
437 if( theNewShape.IsNull() )
440 BOPCol_ListOfShape aShapesList;
441 BOPAlgo_PaveFiller aPaveFiller;
442 HYDROData_MapOfFaceToStricklerType aNewFaces;
444 // add faces to shapes list
445 Iterator anIt( *this );
446 for( ; anIt.More(); anIt.Next() )
447 aShapesList.Append( anIt.Face() );
448 aShapesList.Append( theNewShape );
450 if( aShapesList.Size()==1 && theNewShape.ShapeType()==TopAbs_FACE )
452 aNewFaces.Add( TopoDS::Face( theNewShape ), theNewType );
453 StoreLandCovers( aNewFaces );
457 // prepare pave filler
458 aPaveFiller.SetArguments( aShapesList );
459 aPaveFiller.Perform();
460 Standard_Integer anError = aPaveFiller.ErrorStatus();
464 // add faces to builder
465 BOPAlgo_Builder aBuilder;
467 for( ; anIt.More(); anIt.Next() )
468 aBuilder.AddArgument( anIt.Face() );
469 aBuilder.AddArgument( theNewShape );
471 // perform the partition with the pave filler
472 aBuilder.PerformWithFiller( aPaveFiller );
473 anError = aBuilder.ErrorStatus();
477 // analysis of the history
478 // a. to fill map of shapes which come from the new face
479 NCollection_IndexedMap<TopoDS_Shape> aShapesFromNewFace;
480 //std::cout << "new: " << theNewShape << " " << theNewType << std::endl;
481 TopTools_ListOfShape aModified = aBuilder.Modified( theNewShape );
482 TopTools_ListIteratorOfListOfShape aMIt( aModified );
483 for( ; aMIt.More(); aMIt.Next() )
485 //std::cout << " " << aMIt.Value() << std::endl;
486 aShapesFromNewFace.Add( aMIt.Value() );
489 // b. to fill map of parts except parts from new face
491 for( ; anIt.More(); anIt.Next() )
493 QString aSType = anIt.StricklerType();
494 //std::cout << anIt.Face() << " " << anIt.StricklerType() << std::endl;
495 TopTools_ListOfShape aModified = aBuilder.Modified( anIt.Face() );
496 TopTools_ListIteratorOfListOfShape aMIt( aModified );
497 for( ; aMIt.More(); aMIt.Next() )
499 TopoDS_Shape aShape = aMIt.Value();
500 bool isFace = aShape.ShapeType()==TopAbs_FACE;
501 bool isAlsoFromNew = aShapesFromNewFace.Contains( aShape );
502 //std::cout << " " << aShape << " " << isAlsoFromNew << std::endl;
503 if( isFace && !isAlsoFromNew )
504 aNewFaces.Add( TopoDS::Face( aShape ), aSType );
508 // c. add the new shape if it is face with its type
509 if( theNewShape.ShapeType()==TopAbs_FACE )
510 aNewFaces.Add( TopoDS::Face( theNewShape ), theNewType );
512 // convert map of shape to type to compound and list of types
513 StoreLandCovers( aNewFaces );
518 Replace the set of land covers in the land cover map
519 @param theMap the map of shape (face) to Strickler type (string)
521 void HYDROData_LandCoverMap::StoreLandCovers( const HYDROData_MapOfFaceToStricklerType& theMap )
523 TopoDS_Compound aCompound;
524 BRep_Builder aCompoundBuilder;
525 aCompoundBuilder.MakeCompound( aCompound );
527 int n = theMap.Size();
528 Handle( TDataStd_ExtStringArray ) aTypes =
529 TDataStd_ExtStringArray::Set( myLab.FindChild( DataTag_Types ), 0, n-1, Standard_True );
530 HYDROData_MapOfFaceToStricklerType::Iterator aNFIt( theMap );
531 for( int i=0; aNFIt.More(); aNFIt.Next(), i++ )
533 TopoDS_Face aFace = aNFIt.Key();
534 QString aType = aNFIt.Value();
535 aCompoundBuilder.Add( aCompound, aFace );
536 aTypes->SetValue( i, HYDROData_Tool::toExtString( aType ) );
539 SetShape( aCompound );
543 Find the land cover for the given point
544 @param thePoint the point laying in some land cover
545 @param theType the returned type
546 @return the found land cover's face
548 TopoDS_Face HYDROData_LandCoverMap::FindByPoint( const gp_Pnt2d& thePoint, QString& theType ) const
550 //TODO: some more optimal algorithm
551 Iterator anIt( *this );
552 for( ; anIt.More(); anIt.Next() )
553 if( HYDROData_Tool::ComputePointState( thePoint.XY(), anIt.Face() ) == TopAbs_IN )
555 theType = anIt.StricklerType();
560 return TopoDS_Face();