]> SALOME platform Git repositories - modules/hydro.git/blob - src/HYDROData/HYDROData_ChannelAltitude.cxx
Salome HOME
Merge branch 'BR_H2018_2' of https://codev-tuleap.cea.fr/plugins/git/salome/hydro...
[modules/hydro.git] / src / HYDROData / HYDROData_ChannelAltitude.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_ChannelAltitude.h"
20
21 #include "HYDROData_Document.h"
22 #include "HYDROData_Object.h"
23 #include "HYDROData_Channel.h"
24 #include "HYDROData_Projection.h"
25 #include "HYDROData_Polyline3D.h"
26 #include "HYDROData_PolylineXY.h"
27 #include "HYDROData_ProfileUZ.h"
28 #include "HYDROData_Profile.h"
29
30 //#define _DEVDEBUG_
31 #include "HYDRO_trace.hxx"
32 #include <QString>
33
34 #include <BRepBuilderAPI_MakeVertex.hxx>
35 #include <BRepExtrema_DistShapeShape.hxx>
36
37 #include <Precision.hxx>
38 #include <TopAbs.hxx>
39 #include <TopExp.hxx>
40 #include <TopoDS.hxx>
41 #include <TopoDS_Shape.hxx>
42 #include <TopoDS_Edge.hxx>
43 #include <TopoDS_Face.hxx>
44 #include <TopoDS_Vertex.hxx>
45 #include <TopoDS_Wire.hxx>
46 #include <TopExp_Explorer.hxx>
47 #include <TopAbs_ShapeEnum.hxx>
48 #include <BRep_Tool.hxx>
49 #include <BRepTools.hxx>
50 #include <Geom_Curve.hxx>
51
52 #include <gp_Trsf.hxx>
53 #include <gp_Pnt.hxx>
54 #include <gp_Vec.hxx>
55
56 IMPLEMENT_STANDARD_RTTIEXT(HYDROData_ChannelAltitude, HYDROData_IAltitudeObject)
57
58 HYDROData_ChannelAltitude::HYDROData_ChannelAltitude()
59 : HYDROData_IAltitudeObject()
60 {
61 }
62
63 HYDROData_ChannelAltitude::~HYDROData_ChannelAltitude()
64 {
65 }
66
67 double HYDROData_ChannelAltitude::GetAltitudeForPoint( const gp_XY& thePoint,
68                                                        int theMethod) const
69 {
70   DEBTRACE("GetAltitudeForPoint p(" << thePoint.X() << ", " << thePoint.Y() << ")");
71   double aResAltitude = GetInvalidAltitude();
72
73   Handle(HYDROData_Channel) aChannel =
74   Handle(HYDROData_Channel)::DownCast(GetFatherObject());
75   if (aChannel.IsNull())
76     {
77       DEBTRACE("aChannel.IsNull()");
78       return aResAltitude;
79     }
80   DEBTRACE("aChannel: " << aChannel->GetName().toStdString());
81
82   Handle(HYDROData_Polyline3D) aGuideLine = aChannel->GetGuideLine();
83   if (aGuideLine.IsNull())
84     {
85       DEBTRACE("aGuideLine.IsNull()");
86       return aResAltitude;
87     }
88   //DEBTRACE("aGuideLine: " << aGuideLine->GetName().toStdString());
89
90   Handle(HYDROData_PolylineXY) aGuideXY = aGuideLine->GetPolylineXY();
91   if (aGuideXY.IsNull())
92     {
93       DEBTRACE("aGuideXY.IsNull()");
94       return aResAltitude;
95     }
96   //DEBTRACE("aGuideXY: " << aGuideXY->GetName().toStdString());
97
98   Handle(HYDROData_ProfileUZ) aGuideUZ = aGuideLine->GetProfileUZ();
99   if (aGuideUZ.IsNull())
100     {
101       aGuideUZ = aGuideLine->GetChildProfileUZ(); // profile obtained from bathymetry
102     }
103   if (aGuideUZ.IsNull())
104     {
105       DEBTRACE("aGuideUZ.IsNull()");
106       return aResAltitude;
107     }
108   //DEBTRACE("aGuideUZ: " << aGuideUZ->GetName().toStdString());
109
110   Handle (HYDROData_Profile) aProfile = aChannel->GetProfile();
111   if (aProfile.IsNull())
112     {
113       return aResAltitude;
114     }
115   //DEBTRACE("aProfile: " << aProfile->GetName().toStdString());
116
117   // --- See GEOMImpl_ProjectionDriver.cxx
118
119   TopoDS_Shape aShape =  aGuideXY->GetShape();
120   gp_Pnt P1(thePoint.X(), thePoint.Y(), 0);
121   TopoDS_Shape aPoint = BRepBuilderAPI_MakeVertex(P1).Shape();
122
123   if (aPoint.IsNull() || aShape.IsNull())
124     {
125       DEBTRACE("aPoint.IsNull() || aShape.IsNull()");
126       return aResAltitude;
127     }
128
129   if (aShape.ShapeType() != TopAbs_EDGE && aShape.ShapeType() != TopAbs_WIRE)
130     {
131       DEBTRACE("Projection aborted : the shape is neither an edge nor a wire");
132       return aResAltitude;
133     }
134
135   // Perform projection.
136   BRepExtrema_DistShapeShape aDistShSh(aPoint, aShape, Extrema_ExtFlag_MIN);
137
138   if (aDistShSh.IsDone() == Standard_False)
139     {
140       DEBTRACE("Projection not done");
141       return aResAltitude;
142     }
143
144   Standard_Boolean hasValidSolution = Standard_False;
145   Standard_Integer aNbSolutions = aDistShSh.NbSolution();
146   Standard_Integer i;
147   double aParam = 0.;
148   Standard_Real aTolConf = BRep_Tool::Tolerance(TopoDS::Vertex(aPoint));
149   Standard_Real aTolAng = 1.e-4;
150
151   for (i = 1; i <= aNbSolutions; i++)
152     {
153       Standard_Boolean isValid = Standard_False;
154       BRepExtrema_SupportType aSupportType = aDistShSh.SupportTypeShape2(i);
155       TopoDS_Shape aSupportShape = aDistShSh.SupportOnShape2(i);
156
157       if (aSupportType == BRepExtrema_IsOnEdge)
158         {
159           // Minimal distance inside edge is really a projection.
160           isValid = Standard_True;
161           aDistShSh.ParOnEdgeS2(i, aParam);
162         }
163       else if (aSupportType == BRepExtrema_IsVertex)
164         {
165           TopExp_Explorer anExp(aShape, TopAbs_EDGE);
166
167           if (aDistShSh.Value() <= aTolConf)
168             {
169               // The point lies on the shape. This means this point
170               // is really a projection.
171               for (; anExp.More() && !isValid; anExp.Next())
172                 {
173                   TopoDS_Edge aCurEdge = TopoDS::Edge(anExp.Current());
174
175                   if (aCurEdge.IsNull() == Standard_False)
176                     {
177                       TopoDS_Vertex aVtx[2];
178
179                       TopExp::Vertices(aCurEdge, aVtx[0], aVtx[1]);
180
181                       for (int j = 0; j < 2; j++)
182                         {
183                           if (aSupportShape.IsSame(aVtx[j]))
184                             {
185                               // The current edge is a projection edge.
186                               isValid = Standard_True;
187                               aSupportShape = aCurEdge;
188                               aParam = BRep_Tool::Parameter(aVtx[j], aCurEdge);
189                               break;
190                             }
191                         }
192                     }
193                 }
194             }
195           else
196             {
197               // Minimal distance to vertex is not always a real projection.
198               gp_Pnt aPnt = BRep_Tool::Pnt(TopoDS::Vertex(aPoint));
199               gp_Pnt aPrjPnt = BRep_Tool::Pnt(TopoDS::Vertex(aSupportShape));
200               gp_Vec aDProjP(aPrjPnt, aPnt);
201
202               for (; anExp.More() && !isValid; anExp.Next())
203                 {
204                   TopoDS_Edge aCurEdge = TopoDS::Edge(anExp.Current());
205
206                   if (aCurEdge.IsNull() == Standard_False)
207                     {
208                       TopoDS_Vertex aVtx[2];
209
210                       TopExp::Vertices(aCurEdge, aVtx[0], aVtx[1]);
211
212                       for (int j = 0; j < 2; j++)
213                         {
214                           if (aSupportShape.IsSame(aVtx[j]))
215                             {
216                               // Check if the point is a projection to the current edge.
217                               Standard_Real anEdgePars[2];
218                               Handle(Geom_Curve) aCurve = BRep_Tool::Curve(aCurEdge, anEdgePars[0], anEdgePars[1]);
219                               gp_Pnt aVal;
220                               gp_Vec aD1;
221
222                               aParam = BRep_Tool::Parameter(aVtx[j], aCurEdge);
223                               aCurve->D1(aParam, aVal, aD1);
224
225                               if (Abs(aD1.Dot(aDProjP)) <= aTolAng)
226                                 {
227                                   // The current edge is a projection edge.
228                                   isValid = Standard_True;
229                                   aSupportShape = aCurEdge;
230                                   break;
231                                 }
232                             }
233                         }
234                     }
235                 }
236             }
237         }
238
239       if (isValid)
240         {
241           if (hasValidSolution)
242             {
243               DEBTRACE("Projection aborted : multiple solutions");
244               return aResAltitude;
245             }
246
247           // Store the valid solution.
248           hasValidSolution = Standard_True;
249
250           // profile altitude at projection point
251           HYDROData_IPolyline::PointsList aProfilePoints = aGuideUZ->GetPoints();
252           if ( aProfilePoints.IsEmpty() )
253             {
254               DEBTRACE("empty profile UZ");
255               return aResAltitude;
256             }
257           double aDepth = HYDROData_ProfileUZ::GetDepthFromDistance( aProfilePoints, aParam );
258           //DEBTRACE("profile altitude: " << aDepth);
259
260           // Compute edge index.
261           TopExp_Explorer anExp(aShape, TopAbs_EDGE);
262           int anIndex = 0;
263           for (; anExp.More(); anExp.Next(), anIndex++)
264             {
265               if (aSupportShape.IsSame(anExp.Current()))
266                 {
267                   break;
268                 }
269             }
270
271           // get the XY distance from point to guideline
272           const gp_Pnt &aPntProj = aDistShSh.PointOnShape2(i);
273           //DEBTRACE("projection: (" << aPntProj.X() << ", " << aPntProj.Y() << ", " << aPntProj.Z() << ")");
274           gp_XY aProjXY = gp_XY(aPntProj.X(), aPntProj.Y());
275           aProjXY.Subtract(thePoint);
276           double distance = aProjXY.Modulus();
277           //DEBTRACE("distance to guideline " << distance);
278
279           // get delta altitude on section (supposed symmetric) from guideline distance (aParam)
280           double delta = 0;
281           int i1 = 0;
282           gp_XY pt1 = gp_XY();
283           gp_XY pt2 = gp_XY();
284           HYDROData_ProfileUZ::PointsList aSectionPoints = aProfile->GetParametricPoints();
285           for ( int i = 1, aNbPoints = aSectionPoints.Size(); i <= aNbPoints; ++i )
286             {
287               const HYDROData_IPolyline::Point& aPolylinePoint = aSectionPoints.Value( i );
288               //DEBTRACE("  profile point: " << aPolylinePoint.X() << " " << aPolylinePoint.Y());
289               if (aPolylinePoint.X() < distance)
290                 {
291                   i1 = i;
292                   pt1 = aPolylinePoint;
293                 }
294               if (aPolylinePoint.X() >= distance)
295                 {
296                   pt2 = aPolylinePoint;
297                   break;
298                 }
299             }
300           if ((i1 == 0) && (distance > 0))
301             {
302               DEBTRACE("Projection aborted : non centered profile");
303               return aResAltitude;
304             }
305           if (i1 == aProfilePoints.Size()) // distance >= profile width
306             {
307               delta = pt1.Y();
308             }
309           else
310             {
311               delta = pt1.Y() + (pt2.Y() - pt1.Y())*(distance -pt1.X())/(pt2.X()-pt1.X());
312             }
313           aResAltitude = delta + aDepth;
314           DEBTRACE("distance XY: "<< aParam << " distance to guideline: " << distance << " final altitude: " << aResAltitude << " delta: " << delta);
315           return aResAltitude;
316         }
317     }
318
319   if (!hasValidSolution)
320     {
321       DEBTRACE("Projection aborted : no projection");
322       return aResAltitude;
323     }
324
325   return aResAltitude;
326 }
327
328
329
330