Salome HOME
Dump to python corrected.
[modules/hydro.git] / src / HYDROData / HYDROData_StreamAltitude.cxx
1
2 #include "HYDROData_StreamAltitude.h"
3
4 #include "HYDROData_Document.h"
5 #include "HYDROData_Profile.h"
6 #include "HYDROData_Stream.h"
7 #include "HYDROData_ShapesTool.h"
8
9 #include <BRep_Tool.hxx>
10
11 #include <BRepTopAdaptor_FClass2d.hxx>
12
13 #include <BRepBuilderAPI_MakeEdge.hxx>
14 #include <BRepBuilderAPI_MakeFace.hxx>
15 #include <BRepBuilderAPI_MakeWire.hxx>
16
17 #include <Extrema_ExtElC.hxx>
18
19 #include <GeomAPI_ProjectPointOnCurve.hxx>
20
21 #include <gp_Lin.hxx>
22
23 #include <Precision.hxx>
24
25 #include <TopExp_Explorer.hxx>
26
27 #include <TopoDS.hxx>
28 #include <TopoDS_Wire.hxx>
29 #include <TopoDS_Edge.hxx>
30
31 #include <TopTools_SequenceOfShape.hxx>
32
33 #include <QStringList>
34
35 IMPLEMENT_STANDARD_HANDLE(HYDROData_StreamAltitude, HYDROData_IAltitudeObject)
36 IMPLEMENT_STANDARD_RTTIEXT(HYDROData_StreamAltitude, HYDROData_IAltitudeObject)
37
38 HYDROData_StreamAltitude::HYDROData_StreamAltitude()
39 : HYDROData_IAltitudeObject()
40 {
41 }
42
43 HYDROData_StreamAltitude::~HYDROData_StreamAltitude()
44 {
45 }
46
47 QStringList HYDROData_StreamAltitude::DumpToPython( MapOfTreatedObjects& theTreatedObjects ) const
48 {
49   QStringList aResList = dumpObjectCreation( theTreatedObjects );
50   QString aName = GetObjPyName();
51
52   // TODO
53
54   aResList << QString( "" );
55   aResList << QString( "%1.Update();" ).arg( aName );
56   aResList << QString( "" );
57
58   return aResList;
59 }
60
61 Standard_Real getAltitudeFromProfile( const Handle(HYDROData_Profile)& theProfile,
62                                       const Standard_Real&             theLeftDist,
63                                       const Standard_Real&             theRightDist )
64 {
65   Standard_Real aResAlt = 0.0;
66
67   gp_XY aFirstPoint, aLastPoint;
68   if ( !theProfile->GetLeftPoint( aFirstPoint ) ||
69        !theProfile->GetRightPoint( aLastPoint ) )
70     return aResAlt;
71
72   gp_Pnt aPnt1( aFirstPoint.X(), aFirstPoint.Y(), 0 );
73   gp_Pnt aPnt2( aLastPoint.X(),  aLastPoint.Y(),  0 );
74
75   Standard_Real aProfileDist = aPnt1.Distance( aPnt2 );
76
77   Standard_Real aCoeff = aProfileDist / ( theLeftDist + theRightDist );
78
79   gp_Pnt anIntPoint( aPnt1.XYZ() + ( aCoeff * theLeftDist ) * gp_Dir( gp_Vec( aPnt1, aPnt2 ) ).XYZ() );
80
81   gp_Lin aPointLine( anIntPoint, gp::DZ() );
82
83   gp_Pnt aPrevPoint;
84   gp_Lin aPrevNormal;
85   HYDROData_Profile::ProfilePoints aProfilePoints = theProfile->GetProfilePoints();
86   for ( int i = 1, n = aProfilePoints.Length(); i <= n; ++i )
87   {
88     gp_Pnt aProfPoint( aProfilePoints.Value( i ) );
89
90     Standard_Real aDist = aPointLine.Distance( aProfPoint );
91     if ( aDist <= gp::Resolution() )
92     {
93       // We found the intersected point
94       aResAlt = aProfPoint.Z();
95       break;
96     }
97    
98     gp_Lin aNormal = aPointLine.Normal( aProfPoint );
99     if ( i == 1 )
100     {
101       aPrevNormal = aNormal;
102       aPrevPoint = aProfPoint;
103       continue;
104     }
105
106     if ( aPrevNormal.Direction().Dot( aNormal.Direction() ) < 0 )
107     {
108       // We found the intersected edge
109       gp_Lin anEdgeLine( aPrevPoint, gp_Dir( gp_Vec( aPrevPoint, aProfPoint ) ) );
110      
111       Extrema_ExtElC anExtrema( aPointLine, anEdgeLine, Precision::Angular() );
112       if ( !anExtrema.IsParallel() )
113       {
114         Extrema_POnCurv aFirstPnt, aSecPnt;
115         anExtrema.Points( 1, aFirstPnt, aSecPnt );
116
117         const gp_Pnt& anIntPnt = aSecPnt.Value();
118         aResAlt = anIntPnt.Z();
119
120         break;
121       }
122     }
123
124     aPrevNormal = aNormal;
125     aPrevPoint = aProfPoint;
126   }
127
128   return aResAlt;
129 }
130
131 bool HYDROData_StreamAltitude::getBoundaryProfilesForPoint(
132   const gp_XY&               thePoint,
133   Handle(HYDROData_Profile)& theLeftProfile,
134   Handle(HYDROData_Profile)& theRightProfile ) const
135 {
136   Handle(HYDROData_Stream) aStream =
137     Handle(HYDROData_Stream)::DownCast( GetFatherObject() );
138   if ( aStream.IsNull() )
139     return false;
140
141   HYDROData_SequenceOfObjects aStreamProfiles = aStream->GetProfiles();
142   if ( aStreamProfiles.Length() < 2 )
143     return false;
144
145   Handle(HYDROData_Profile) aPrevProfile;
146   gp_Pnt aPrevPnt1, aPrevPnt2;
147   for ( int i = 1, n = aStreamProfiles.Length(); i <= n; ++i )
148   {
149     Handle(HYDROData_Profile) aProfile =
150       Handle(HYDROData_Profile)::DownCast( aStreamProfiles.Value( i ) );
151     if ( aProfile.IsNull() )
152       continue;
153
154     gp_XY aFirstPoint, aLastPoint;
155     if ( !aProfile->GetLeftPoint( aFirstPoint ) || !aProfile->GetRightPoint( aLastPoint ) )
156       continue;
157
158     gp_Pnt aPnt1( aFirstPoint.X(), aFirstPoint.Y(), 0 );
159     gp_Pnt aPnt2( aLastPoint.X(),  aLastPoint.Y(),  0 );
160
161     if ( !aPrevProfile.IsNull() )
162     {
163       BRepBuilderAPI_MakeEdge aLeftMakeEdge( aPrevPnt1, aPrevPnt2 );
164       BRepBuilderAPI_MakeEdge aBotMakeEdge( aPrevPnt2, aPnt2 );
165       BRepBuilderAPI_MakeEdge aRightMakeEdge( aPnt2, aPnt1 );
166       BRepBuilderAPI_MakeEdge aTopMakeEdge( aPnt1, aPrevPnt1 );
167
168       BRepBuilderAPI_MakeWire aMakeWire( aLeftMakeEdge.Edge(), aBotMakeEdge.Edge(), 
169                                          aRightMakeEdge.Edge(), aTopMakeEdge.Edge() );
170
171       BRepBuilderAPI_MakeFace aMakeFace( aMakeWire.Wire() );
172
173       TopoDS_Face aProfilesFace = aMakeFace.Face();
174
175       BRepTopAdaptor_FClass2d aClassifier( aProfilesFace, Precision::Confusion() );
176       TopAbs_State aPointState = aClassifier.Perform( gp_Pnt2d( thePoint ), Standard_False );
177       if ( aPointState != TopAbs_OUT )
178       {
179         theLeftProfile = aPrevProfile;
180         theRightProfile = aProfile;
181         break;
182       }
183     }
184
185     aPrevProfile = aProfile;
186     aPrevPnt1 = aPnt1;
187     aPrevPnt2 = aPnt2;
188   }
189
190   return !theLeftProfile.IsNull() && !theRightProfile.IsNull();
191 }
192
193 double HYDROData_StreamAltitude::GetAltitudeForPoint( const gp_XY& thePoint ) const
194 {
195   double aResAltitude = GetInvalidAltitude();
196
197   Handle(HYDROData_Stream) aStream =
198     Handle(HYDROData_Stream)::DownCast( GetFatherObject() );
199   if ( aStream.IsNull() )
200     return aResAltitude;
201
202   TopoDS_Shape aStreamShape = aStream->GetTopShape();
203   if ( aStreamShape.IsNull() )
204     return aResAltitude;
205
206   TopExp_Explorer aStreamFaceExp( aStreamShape, TopAbs_FACE );
207   if ( !aStreamFaceExp.More() )
208     return aResAltitude;
209
210   // Get only face because of 2d profile wires is in compound
211   TopoDS_Face aStreamFace = TopoDS::Face( aStreamFaceExp.Current() );
212
213   // Check if point is inside of stream presentation
214   BRepTopAdaptor_FClass2d aClassifier( aStreamFace, Precision::Confusion() );
215   TopAbs_State aPointState = aClassifier.Perform( gp_Pnt2d( thePoint ), Standard_False );
216   if ( aPointState == TopAbs_OUT )
217     return aResAltitude;
218
219   // Find the two profiles between which the point is lies
220   Handle(HYDROData_Profile) aLeftProfile, aRightProfile;
221   if ( !getBoundaryProfilesForPoint( thePoint, aLeftProfile, aRightProfile ) )
222     return aResAltitude;
223
224   // Find the projections of point to borders of stream
225   gp_XYZ aPointToTest( thePoint.X(), thePoint.Y(), 0.0 );
226
227   TopoDS_Edge aStreamLeftEdge = TopoDS::Edge( aStream->GetLeftShape() );
228   TopoDS_Edge aStreamRightEdge = TopoDS::Edge( aStream->GetRightShape() );
229
230   Standard_Real aFirst = 0.0, aLast = 0.0;
231   Handle(Geom_Curve) anEdgeLeftCurve = BRep_Tool::Curve( aStreamLeftEdge, aFirst, aLast );
232   Handle(Geom_Curve) anEdgeRightCurve = BRep_Tool::Curve( aStreamRightEdge, aFirst, aLast );
233
234   GeomAPI_ProjectPointOnCurve aLeftProject( aPointToTest, anEdgeLeftCurve );
235   GeomAPI_ProjectPointOnCurve aRightProject( aPointToTest, anEdgeRightCurve );
236
237   Standard_Real aLeftDist = aLeftProject.LowerDistance();
238   Standard_Real aRightDist = aRightProject.LowerDistance();
239
240   // Find the altitude in profiles
241   Standard_Real aLeftAlt = getAltitudeFromProfile( aLeftProfile, aLeftDist, aRightDist );
242   Standard_Real aRightAlt = getAltitudeFromProfile( aRightProfile, aLeftDist, aRightDist );
243
244   // Interpolate altitudes
245   Standard_Real aFirstCoeff = aLeftDist / ( aLeftDist + aRightDist );
246   Standard_Real aSecCoeff = aRightDist / ( aLeftDist + aRightDist );
247
248   aResAltitude = aLeftAlt * aFirstCoeff + aRightAlt * aSecCoeff;
249
250   return aResAltitude;
251 }
252
253
254
255