#include "HYDROData_StreamAltitude.h"
+
#include "HYDROData_Document.h"
+#include "HYDROData_Profile.h"
+#include "HYDROData_Stream.h"
+#include "HYDROData_ShapesTool.h"
+
+#include <BRep_Tool.hxx>
+
+#include <BRepTopAdaptor_FClass2d.hxx>
+
+#include <BRepBuilderAPI_MakeEdge.hxx>
+#include <BRepBuilderAPI_MakeFace.hxx>
+#include <BRepBuilderAPI_MakeWire.hxx>
+
+#include <Extrema_ExtElC.hxx>
+
+#include <GeomAPI_ProjectPointOnCurve.hxx>
+
+#include <gp_Lin.hxx>
+
+#include <Precision.hxx>
+
+#include <TopExp_Explorer.hxx>
+
+#include <TopoDS.hxx>
+#include <TopoDS_Wire.hxx>
+#include <TopoDS_Edge.hxx>
+
+#include <TopTools_SequenceOfShape.hxx>
#include <QStringList>
return aResList;
}
+Standard_Real getAltitudeFromProfile( const Handle(HYDROData_Profile)& theProfile,
+ const Standard_Real& theLeftDist,
+ const Standard_Real& theRightDist )
+{
+ Standard_Real aResAlt = 0.0;
+
+ gp_XY aFirstPoint, aLastPoint;
+ if ( !theProfile->GetLeftPoint( aFirstPoint ) ||
+ !theProfile->GetRightPoint( aLastPoint ) )
+ return aResAlt;
+
+ gp_Pnt aPnt1( aFirstPoint.X(), aFirstPoint.Y(), 0 );
+ gp_Pnt aPnt2( aLastPoint.X(), aLastPoint.Y(), 0 );
+
+ Standard_Real aProfileDist = aPnt1.Distance( aPnt2 );
+
+ Standard_Real aCoeff = aProfileDist / ( theLeftDist + theRightDist );
+
+ gp_Pnt anIntPoint( aPnt1.XYZ() + ( aCoeff * theLeftDist ) * gp_Dir( gp_Vec( aPnt1, aPnt2 ) ).XYZ() );
+
+ gp_Lin aPointLine( anIntPoint, gp::DZ() );
+
+ gp_Pnt aPrevPoint;
+ gp_Lin aPrevNormal;
+ HYDROData_Profile::ProfilePoints aProfilePoints = theProfile->GetProfilePoints();
+ for ( int i = 1, n = aProfilePoints.Length(); i <= n; ++i )
+ {
+ gp_Pnt aProfPoint( aProfilePoints.Value( i ) );
+
+ Standard_Real aDist = aPointLine.Distance( aProfPoint );
+ if ( aDist <= gp::Resolution() )
+ {
+ // We found the intersected point
+ aResAlt = aProfPoint.Z();
+ break;
+ }
+
+ gp_Lin aNormal = aPointLine.Normal( aProfPoint );
+ if ( i == 1 )
+ {
+ aPrevNormal = aNormal;
+ aPrevPoint = aProfPoint;
+ continue;
+ }
+
+ if ( aPrevNormal.Direction().Dot( aNormal.Direction() ) < 0 )
+ {
+ // We found the intersected edge
+ gp_Lin anEdgeLine( aPrevPoint, gp_Dir( gp_Vec( aPrevPoint, aProfPoint ) ) );
+
+ Extrema_ExtElC anExtrema( aPointLine, anEdgeLine, Precision::Angular() );
+ if ( !anExtrema.IsParallel() )
+ {
+ Extrema_POnCurv aFirstPnt, aSecPnt;
+ anExtrema.Points( 1, aFirstPnt, aSecPnt );
+
+ const gp_Pnt& anIntPnt = aSecPnt.Value();
+ aResAlt = anIntPnt.Z();
+
+ break;
+ }
+ }
+
+ aPrevNormal = aNormal;
+ aPrevPoint = aProfPoint;
+ }
+
+ return aResAlt;
+}
+
+bool HYDROData_StreamAltitude::getBoundaryProfilesForPoint(
+ const gp_XY& thePoint,
+ Handle(HYDROData_Profile)& theLeftProfile,
+ Handle(HYDROData_Profile)& theRightProfile ) const
+{
+ Handle(HYDROData_Stream) aStream =
+ Handle(HYDROData_Stream)::DownCast( GetFatherObject() );
+ if ( aStream.IsNull() )
+ return false;
+
+ HYDROData_SequenceOfObjects aStreamProfiles = aStream->GetProfiles();
+ if ( aStreamProfiles.Length() < 2 )
+ return false;
+
+ Handle(HYDROData_Profile) aPrevProfile;
+ gp_Pnt aPrevPnt1, aPrevPnt2;
+ for ( int i = 1, n = aStreamProfiles.Length(); i <= n; ++i )
+ {
+ Handle(HYDROData_Profile) aProfile =
+ Handle(HYDROData_Profile)::DownCast( aStreamProfiles.Value( i ) );
+ if ( aProfile.IsNull() )
+ continue;
+
+ gp_XY aFirstPoint, aLastPoint;
+ if ( !aProfile->GetLeftPoint( aFirstPoint ) || !aProfile->GetRightPoint( aLastPoint ) )
+ continue;
+
+ gp_Pnt aPnt1( aFirstPoint.X(), aFirstPoint.Y(), 0 );
+ gp_Pnt aPnt2( aLastPoint.X(), aLastPoint.Y(), 0 );
+
+ if ( !aPrevProfile.IsNull() )
+ {
+ BRepBuilderAPI_MakeEdge aLeftMakeEdge( aPrevPnt1, aPrevPnt2 );
+ BRepBuilderAPI_MakeEdge aBotMakeEdge( aPrevPnt2, aPnt2 );
+ BRepBuilderAPI_MakeEdge aRightMakeEdge( aPnt2, aPnt1 );
+ BRepBuilderAPI_MakeEdge aTopMakeEdge( aPnt1, aPrevPnt1 );
+
+ BRepBuilderAPI_MakeWire aMakeWire( aLeftMakeEdge.Edge(), aBotMakeEdge.Edge(),
+ aRightMakeEdge.Edge(), aTopMakeEdge.Edge() );
+
+ BRepBuilderAPI_MakeFace aMakeFace( aMakeWire.Wire() );
+
+ TopoDS_Face aProfilesFace = aMakeFace.Face();
+
+ BRepTopAdaptor_FClass2d aClassifier( aProfilesFace, Precision::Confusion() );
+ TopAbs_State aPointState = aClassifier.Perform( gp_Pnt2d( thePoint ), Standard_False );
+ if ( aPointState != TopAbs_OUT )
+ {
+ theLeftProfile = aPrevProfile;
+ theRightProfile = aProfile;
+ break;
+ }
+ }
+
+ aPrevProfile = aProfile;
+ aPrevPnt1 = aPnt1;
+ aPrevPnt2 = aPnt2;
+ }
+
+ return !theLeftProfile.IsNull() && !theRightProfile.IsNull();
+}
+
double HYDROData_StreamAltitude::GetAltitudeForPoint( const gp_XY& thePoint ) const
{
double aResAltitude = GetInvalidAltitude();
+ Handle(HYDROData_Stream) aStream =
+ Handle(HYDROData_Stream)::DownCast( GetFatherObject() );
+ if ( aStream.IsNull() )
+ return aResAltitude;
+
+ TopoDS_Shape aStreamShape = aStream->GetTopShape();
+ if ( aStreamShape.IsNull() )
+ return aResAltitude;
+
+ TopExp_Explorer aStreamFaceExp( aStreamShape, TopAbs_FACE );
+ if ( !aStreamFaceExp.More() )
+ return aResAltitude;
+
+ // Get only face because of 2d profile wires is in compound
+ TopoDS_Face aStreamFace = TopoDS::Face( aStreamFaceExp.Current() );
+
+ // Check if point is inside of stream presentation
+ BRepTopAdaptor_FClass2d aClassifier( aStreamFace, Precision::Confusion() );
+ TopAbs_State aPointState = aClassifier.Perform( gp_Pnt2d( thePoint ), Standard_False );
+ if ( aPointState == TopAbs_OUT )
+ return aResAltitude;
+
+ // Find the two profiles between which the point is lies
+ Handle(HYDROData_Profile) aLeftProfile, aRightProfile;
+ if ( !getBoundaryProfilesForPoint( thePoint, aLeftProfile, aRightProfile ) )
+ return aResAltitude;
+
+ // Find the projections of point to borders of stream
+ gp_XYZ aPointToTest( thePoint.X(), thePoint.Y(), 0.0 );
+
+ TopoDS_Edge aStreamLeftEdge = TopoDS::Edge( aStream->GetLeftShape() );
+ TopoDS_Edge aStreamRightEdge = TopoDS::Edge( aStream->GetRightShape() );
+
+ Standard_Real aFirst = 0.0, aLast = 0.0;
+ Handle(Geom_Curve) anEdgeLeftCurve = BRep_Tool::Curve( aStreamLeftEdge, aFirst, aLast );
+ Handle(Geom_Curve) anEdgeRightCurve = BRep_Tool::Curve( aStreamRightEdge, aFirst, aLast );
+
+ GeomAPI_ProjectPointOnCurve aLeftProject( aPointToTest, anEdgeLeftCurve );
+ GeomAPI_ProjectPointOnCurve aRightProject( aPointToTest, anEdgeRightCurve );
+
+ Standard_Real aLeftDist = aLeftProject.LowerDistance();
+ Standard_Real aRightDist = aRightProject.LowerDistance();
+
+ // Find the altitude in profiles
+ Standard_Real aLeftAlt = getAltitudeFromProfile( aLeftProfile, aLeftDist, aRightDist );
+ Standard_Real aRightAlt = getAltitudeFromProfile( aRightProfile, aLeftDist, aRightDist );
+
+ // Interpolate altitudes
+ Standard_Real aFirstCoeff = aLeftDist / ( aLeftDist + aRightDist );
+ Standard_Real aSecCoeff = aRightDist / ( aLeftDist + aRightDist );
+
+ aResAltitude = aLeftAlt * aFirstCoeff + aRightAlt * aSecCoeff;
+
return aResAltitude;
}
