1 // Copyright (C) 2007-2015 CEA/DEN, EDF R&D, OPEN CASCADE
3 // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
23 #include "HYDROData_StreamAltitude.h"
25 #include "HYDROData_Document.h"
26 #include "HYDROData_Profile.h"
27 #include "HYDROData_Stream.h"
28 #include "HYDROData_ShapesTool.h"
29 #include <HYDROData_Tool.h>
30 #include <BRep_Tool.hxx>
32 #include <BRepBuilderAPI_MakeEdge.hxx>
33 #include <BRepBuilderAPI_MakeFace.hxx>
34 #include <BRepBuilderAPI_MakeWire.hxx>
36 #include <Extrema_ExtElC.hxx>
38 #include <GeomAPI_ProjectPointOnCurve.hxx>
42 #include <Precision.hxx>
44 #include <TopExp_Explorer.hxx>
47 #include <TopoDS_Wire.hxx>
48 #include <TopoDS_Edge.hxx>
50 #include <TopTools_SequenceOfShape.hxx>
52 #include <Geom_Line.hxx>
54 #include <QStringList>
57 #include <BRepTools.hxx>
58 #include <BRep_Builder.hxx>
59 #include <BRepBuilderAPI_MakeVertex.hxx>
61 IMPLEMENT_STANDARD_HANDLE(HYDROData_StreamAltitude, HYDROData_IAltitudeObject)
62 IMPLEMENT_STANDARD_RTTIEXT(HYDROData_StreamAltitude, HYDROData_IAltitudeObject)
64 HYDROData_StreamAltitude::HYDROData_StreamAltitude()
65 : HYDROData_IAltitudeObject()
69 HYDROData_StreamAltitude::~HYDROData_StreamAltitude()
73 Standard_Real getAltitudeFromProfile( const Handle(HYDROData_Profile)& theProfile,
74 const Standard_Real& theLeftDist,
75 const Standard_Real& theRightDist )
77 Standard_Real aResAlt = 0.0;
79 gp_XY aFirstPoint, aLastPoint;
80 if ( !theProfile->GetLeftPoint( aFirstPoint, false ) ||
81 !theProfile->GetRightPoint( aLastPoint, false ) )
84 gp_Pnt aPnt1( aFirstPoint.X(), aFirstPoint.Y(), 0 );
85 gp_Pnt aPnt2( aLastPoint.X(), aLastPoint.Y(), 0 );
87 Standard_Real aProfileDist = aPnt1.Distance( aPnt2 );
89 Standard_Real aCoeff = aProfileDist / ( theLeftDist + theRightDist );
91 gp_Pnt anIntPoint( aPnt1.XYZ() + ( aCoeff * theLeftDist ) * gp_Dir( gp_Vec( aPnt1, aPnt2 ) ).XYZ() );
93 gp_Lin aPointLine( anIntPoint, gp::DZ() );
97 HYDROData_Profile::ProfilePoints aProfilePoints = theProfile->GetProfilePoints( false );
98 for ( int i = 1, n = aProfilePoints.Length(); i <= n; ++i )
100 gp_Pnt aProfPoint( aProfilePoints.Value( i ) );
102 Standard_Real aDist = aPointLine.Distance( aProfPoint );
103 if ( aDist <= gp::Resolution() )
105 // We found the intersected point
106 aResAlt = aProfPoint.Z();
110 gp_Lin aNormal = aPointLine.Normal( aProfPoint );
113 aPrevNormal = aNormal;
114 aPrevPoint = aProfPoint;
118 if ( aPrevNormal.Direction().Dot( aNormal.Direction() ) < 0 )
120 // We found the intersected edge
121 gp_Lin anEdgeLine( aPrevPoint, gp_Dir( gp_Vec( aPrevPoint, aProfPoint ) ) );
123 Extrema_ExtElC anExtrema( aPointLine, anEdgeLine, Precision::Angular() );
124 if ( !anExtrema.IsParallel() )
126 Extrema_POnCurv aFirstPnt, aSecPnt;
127 anExtrema.Points( 1, aFirstPnt, aSecPnt );
129 const gp_Pnt& anIntPnt = aSecPnt.Value();
130 aResAlt = anIntPnt.Z();
136 aPrevNormal = aNormal;
137 aPrevPoint = aProfPoint;
143 bool HYDROData_StreamAltitude::getBoundaryProfilesForPoint(
144 const gp_XY& thePoint,
145 Handle(HYDROData_Profile)& theLeftProfile,
146 Handle(HYDROData_Profile)& theRightProfile ) const
148 Handle(HYDROData_Stream) aStream =
149 Handle(HYDROData_Stream)::DownCast( GetFatherObject() );
150 if ( aStream.IsNull() )
153 HYDROData_SequenceOfObjects aStreamProfiles = aStream->GetProfiles();
154 if ( aStreamProfiles.Length() < 2 )
157 Handle(HYDROData_Profile) aPrevProfile;
158 gp_Pnt aPrevPnt1, aPrevPnt2;
159 for ( int i = 1, n = aStreamProfiles.Length(); i <= n; ++i )
161 Handle(HYDROData_Profile) aProfile =
162 Handle(HYDROData_Profile)::DownCast( aStreamProfiles.Value( i ) );
163 if ( aProfile.IsNull() )
166 gp_XY aFirstPoint, aLastPoint;
167 if ( !aProfile->GetLeftPoint( aFirstPoint, false ) ||
168 !aProfile->GetRightPoint( aLastPoint, false ) )
171 gp_Pnt aPnt1( aFirstPoint.X(), aFirstPoint.Y(), 0 );
172 gp_Pnt aPnt2( aLastPoint.X(), aLastPoint.Y(), 0 );
174 if ( !aPrevProfile.IsNull() )
176 BRepBuilderAPI_MakeEdge aLeftMakeEdge( aPrevPnt1, aPrevPnt2 );
177 BRepBuilderAPI_MakeEdge aBotMakeEdge( aPrevPnt2, aPnt2 );
178 BRepBuilderAPI_MakeEdge aRightMakeEdge( aPnt2, aPnt1 );
179 BRepBuilderAPI_MakeEdge aTopMakeEdge( aPnt1, aPrevPnt1 );
181 BRepBuilderAPI_MakeWire aMakeWire( aLeftMakeEdge.Edge(), aBotMakeEdge.Edge(),
182 aRightMakeEdge.Edge(), aTopMakeEdge.Edge() );
184 BRepBuilderAPI_MakeFace aMakeFace( aMakeWire.Wire() );
186 TopoDS_Face aProfilesFace = aMakeFace.Face();
188 TopoDS_Compound aCmp;
190 aBB.MakeCompound(aCmp);
191 aBB.Add(aCmp, aProfilesFace);
192 gp_Pnt aPnt (thePoint.X(), thePoint.Y(), 0.);
193 BRepBuilderAPI_MakeVertex aMk(aPnt);
194 aBB.Add(aCmp, aMk.Vertex());
195 BRepTools::Write(aCmp, "ProfileFace.brep");
198 TopAbs_State aPointState = HYDROData_Tool::ComputePointState(thePoint, aProfilesFace);
201 cout << "Point status is = " << aPointState <<endl;
203 if ( aPointState != TopAbs_OUT )
205 theLeftProfile = aPrevProfile;
206 theRightProfile = aProfile;
211 aPrevProfile = aProfile;
216 return !theLeftProfile.IsNull() && !theRightProfile.IsNull();
219 double HYDROData_StreamAltitude::GetAltitudeForPoint( const gp_XY& thePoint ) const
221 double aResAltitude = GetInvalidAltitude();
223 Handle(HYDROData_Stream) aStream =
224 Handle(HYDROData_Stream)::DownCast( GetFatherObject() );
225 if ( aStream.IsNull() )
228 TopoDS_Shape aStreamShape = aStream->GetTopShape();
229 if ( aStreamShape.IsNull() )
232 TopExp_Explorer aStreamFaceExp( aStreamShape, TopAbs_FACE );
233 if ( !aStreamFaceExp.More() )
236 // Get only face because of 2d profile wires is in compound
237 TopoDS_Face aStreamFace = TopoDS::Face( aStreamFaceExp.Current() );
239 // Check if point is inside of stream presentation
240 TopAbs_State aPointState = HYDROData_Tool::ComputePointState(thePoint, aStreamFace);
243 cout << "Point status is = " << aPointState <<endl;
244 TopoDS_Compound aCmp;
246 aBB.MakeCompound(aCmp);
247 aBB.Add(aCmp, aStreamFace);
248 gp_Pnt aPnt (thePoint.X(), thePoint.Y(), 0.);
249 BRepBuilderAPI_MakeVertex aMk(aPnt);
250 aBB.Add(aCmp, aMk.Vertex());
251 BRepTools::Write(aCmp, "FCL2d.brep");
253 if ( aPointState == TopAbs_OUT )
256 // Find the two profiles between which the point is lies
257 Handle(HYDROData_Profile) aLeftProfile, aRightProfile;
258 if ( !getBoundaryProfilesForPoint( thePoint, aLeftProfile, aRightProfile ) )
261 // Find the projections of point to borders of stream
262 gp_XYZ aPointToTest( thePoint.X(), thePoint.Y(), 0.0 );
264 TopoDS_Edge aStreamLeftEdge = TopoDS::Edge( aStream->GetLeftShape() );
265 TopoDS_Edge aStreamRightEdge = TopoDS::Edge( aStream->GetRightShape() );
267 Standard_Real aFirst = 0.0, aLast = 0.0;
268 Handle(Geom_Curve) anEdgeLeftCurve = BRep_Tool::Curve( aStreamLeftEdge, aFirst, aLast );
269 Handle(Geom_Curve) anEdgeRightCurve = BRep_Tool::Curve( aStreamRightEdge, aFirst, aLast );
271 GeomAPI_ProjectPointOnCurve aLeftProject( aPointToTest, anEdgeLeftCurve );
272 GeomAPI_ProjectPointOnCurve aRightProject( aPointToTest, anEdgeRightCurve );
274 Standard_Real aLeftDist = aLeftProject.LowerDistance();
275 Standard_Real aRightDist = aRightProject.LowerDistance();
277 // Find the altitude in profiles
278 Standard_Real aLeftAlt = getAltitudeFromProfile( aLeftProfile, aLeftDist, aRightDist );
279 Standard_Real aRightAlt = getAltitudeFromProfile( aRightProfile, aLeftDist, aRightDist );
281 // Interpolate altitudes
282 // Left profile line ( the segment between the firts and the last profile point )
283 HYDROData_Profile::ProfilePoints aLeftProfilePoints = aLeftProfile->GetProfilePoints( false );
284 gp_Pnt aLeftProfileP1( aLeftProfilePoints.First() );
285 aLeftProfileP1.SetZ( 0 );
286 gp_Pnt aLeftProfileP2( aLeftProfilePoints.Last() );
287 aLeftProfileP2.SetZ( 0 );
288 gp_Vec aLeftProfileVec( aLeftProfileP1, aLeftProfileP2 );
289 Handle(Geom_Line) aLeftProfileLine = new Geom_Line( gp_Ax1( aLeftProfileP1, aLeftProfileVec ) );
290 // Right profile line
291 HYDROData_Profile::ProfilePoints aRightProfilePoints = aRightProfile->GetProfilePoints( false );
292 gp_Pnt aRightProfileP1( aRightProfilePoints.First() );
293 aRightProfileP1.SetZ( 0 );
294 gp_Pnt aRightProfileP2( aRightProfilePoints.Last() );
295 aRightProfileP2.SetZ( 0 );
296 gp_Vec aRightProfileVec( aRightProfileP1, aRightProfileP2 );
297 Handle(Geom_Line) aRightProfileLine = new Geom_Line( gp_Ax1( aRightProfileP1, aRightProfileVec ) );
298 // The point projections on the left and right profiles
299 GeomAPI_ProjectPointOnCurve aLeftProfileProject( aPointToTest, aLeftProfileLine );
300 GeomAPI_ProjectPointOnCurve aRightProfileProject( aPointToTest, aRightProfileLine );
301 // The point distance to the left and right profiles
302 Standard_Real aLeftProfileDist = aLeftProfileProject.LowerDistance();
303 Standard_Real aRightProfileDist = aRightProfileProject.LowerDistance();
305 Standard_Real aFirstCoeff = aLeftProfileDist / ( aLeftProfileDist + aRightProfileDist );
306 Standard_Real aSecCoeff = aRightProfileDist / ( aLeftProfileDist + aRightProfileDist );
308 aResAltitude = aLeftAlt * ( 1 - aFirstCoeff ) + aRightAlt * ( 1 - aSecCoeff );