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_ChannelAltitude.h"
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"
31 #include "HYDRO_trace.hxx"
34 #include <BRepBuilderAPI_MakeVertex.hxx>
35 #include <BRepExtrema_DistShapeShape.hxx>
37 #include <Precision.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>
52 #include <gp_Trsf.hxx>
56 IMPLEMENT_STANDARD_HANDLE(HYDROData_ChannelAltitude, HYDROData_IAltitudeObject)
57 IMPLEMENT_STANDARD_RTTIEXT(HYDROData_ChannelAltitude, HYDROData_IAltitudeObject)
59 HYDROData_ChannelAltitude::HYDROData_ChannelAltitude()
60 : HYDROData_IAltitudeObject()
64 HYDROData_ChannelAltitude::~HYDROData_ChannelAltitude()
68 double HYDROData_ChannelAltitude::GetAltitudeForPoint( const gp_XY& thePoint ) const
70 DEBTRACE("HYDROData_ChannelAltitude::GetAltitudeForPoint");
71 double aResAltitude = GetInvalidAltitude();
73 Handle(HYDROData_Channel) aChannel =
74 Handle(HYDROData_Channel)::DownCast(GetFatherObject());
75 if (aChannel.IsNull())
77 DEBTRACE("aChannel.IsNull()");
80 DEBTRACE("aChannel: " << aChannel->GetName().toStdString());
82 Handle(HYDROData_Polyline3D) aGuideLine = aChannel->GetGuideLine();
83 if (aGuideLine.IsNull())
85 DEBTRACE("aGuideLine.IsNull()");
88 DEBTRACE("aGuideLine: " << aGuideLine->GetName().toStdString());
90 Handle(HYDROData_PolylineXY) aGuideXY = aGuideLine->GetPolylineXY();
91 if (aGuideXY.IsNull())
93 DEBTRACE("aGuideXY.IsNull()");
96 DEBTRACE("aGuideXY: " << aGuideXY->GetName().toStdString());
98 Handle(HYDROData_ProfileUZ) aGuideUZ = aGuideLine->GetProfileUZ();
99 if (aGuideUZ.IsNull())
101 aGuideUZ = aGuideLine->GetChildProfileUZ(); // profile obtained from bathymetry
103 if (aGuideUZ.IsNull())
105 DEBTRACE("aGuideUZ.IsNull()");
108 DEBTRACE("aGuideUZ: " << aGuideUZ->GetName().toStdString());
110 Handle (HYDROData_Profile) aProfile = aChannel->GetProfile();
111 if (aProfile.IsNull())
115 DEBTRACE("aProfile: " << aProfile->GetName().toStdString());
117 // --- See GEOMImpl_ProjectionDriver.cxx
119 // TopoDS_Shape aShape = aGuideLine->GetShape3D();
120 // double middleZ = -9999;
121 // aGuideLine->GetMiddleZ(middleZ); // use the middle Z value of the 3d line to help the projection.
122 // if (middleZ < -9000)
124 // DEBTRACE("the middle Z value of the 3d line is incorrect");
126 // gp_Pnt P1(thePoint.X(), thePoint.Y(), middleZ);
128 TopoDS_Shape aShape = aGuideXY->GetShape();
129 gp_Pnt P1(thePoint.X(), thePoint.Y(), 0);
130 TopoDS_Shape aPoint = BRepBuilderAPI_MakeVertex(P1).Shape();
132 if (aPoint.IsNull() || aShape.IsNull())
134 DEBTRACE("aPoint.IsNull() || aShape.IsNull()");
138 if (aShape.ShapeType() != TopAbs_EDGE && aShape.ShapeType() != TopAbs_WIRE)
140 DEBTRACE("Projection aborted : the shape is neither an edge nor a wire");
144 // Perform projection.
145 BRepExtrema_DistShapeShape aDistShSh(aPoint, aShape, Extrema_ExtFlag_MIN);
147 if (aDistShSh.IsDone() == Standard_False)
149 DEBTRACE("Projection not done");
153 Standard_Boolean hasValidSolution = Standard_False;
154 Standard_Integer aNbSolutions = aDistShSh.NbSolution();
157 Standard_Real aTolConf = BRep_Tool::Tolerance(TopoDS::Vertex(aPoint));
158 Standard_Real aTolAng = 1.e-4;
160 for (i = 1; i <= aNbSolutions; i++)
162 Standard_Boolean isValid = Standard_False;
163 BRepExtrema_SupportType aSupportType = aDistShSh.SupportTypeShape2(i);
164 TopoDS_Shape aSupportShape = aDistShSh.SupportOnShape2(i);
166 if (aSupportType == BRepExtrema_IsOnEdge)
168 // Minimal distance inside edge is really a projection.
169 isValid = Standard_True;
170 aDistShSh.ParOnEdgeS2(i, aParam);
172 else if (aSupportType == BRepExtrema_IsVertex)
174 TopExp_Explorer anExp(aShape, TopAbs_EDGE);
176 if (aDistShSh.Value() <= aTolConf)
178 // The point lies on the shape. This means this point
179 // is really a projection.
180 for (; anExp.More() && !isValid; anExp.Next())
182 TopoDS_Edge aCurEdge = TopoDS::Edge(anExp.Current());
184 if (aCurEdge.IsNull() == Standard_False)
186 TopoDS_Vertex aVtx[2];
188 TopExp::Vertices(aCurEdge, aVtx[0], aVtx[1]);
190 for (int j = 0; j < 2; j++)
192 if (aSupportShape.IsSame(aVtx[j]))
194 // The current edge is a projection edge.
195 isValid = Standard_True;
196 aSupportShape = aCurEdge;
197 aParam = BRep_Tool::Parameter(aVtx[j], aCurEdge);
206 // Minimal distance to vertex is not always a real projection.
207 gp_Pnt aPnt = BRep_Tool::Pnt(TopoDS::Vertex(aPoint));
208 gp_Pnt aPrjPnt = BRep_Tool::Pnt(TopoDS::Vertex(aSupportShape));
209 gp_Vec aDProjP(aPrjPnt, aPnt);
211 for (; anExp.More() && !isValid; anExp.Next())
213 TopoDS_Edge aCurEdge = TopoDS::Edge(anExp.Current());
215 if (aCurEdge.IsNull() == Standard_False)
217 TopoDS_Vertex aVtx[2];
219 TopExp::Vertices(aCurEdge, aVtx[0], aVtx[1]);
221 for (int j = 0; j < 2; j++)
223 if (aSupportShape.IsSame(aVtx[j]))
225 // Check if the point is a projection to the current edge.
226 Standard_Real anEdgePars[2];
227 Handle(Geom_Curve) aCurve = BRep_Tool::Curve(aCurEdge, anEdgePars[0], anEdgePars[1]);
231 aParam = BRep_Tool::Parameter(aVtx[j], aCurEdge);
232 aCurve->D1(aParam, aVal, aD1);
234 if (Abs(aD1.Dot(aDProjP)) <= aTolAng)
236 // The current edge is a projection edge.
237 isValid = Standard_True;
238 aSupportShape = aCurEdge;
250 if (hasValidSolution)
252 DEBTRACE("Projection aborted : multiple solutions");
256 // Store the valid solution.
257 hasValidSolution = Standard_True;
259 // Normalize parameter.
260 TopoDS_Edge aSupportEdge = TopoDS::Edge(aSupportShape);
261 Standard_Real aF, aL;
262 BRep_Tool::Range(aSupportEdge, aF, aL);
263 if (Abs(aL - aF) <= aTolConf)
265 DEBTRACE("Projection aborted : degenerated projection edge");
268 aParam = (aParam - aF)/(aL - aF);
270 // Compute edge index.
271 TopExp_Explorer anExp(aShape, TopAbs_EDGE);
273 for (; anExp.More(); anExp.Next(), anIndex++)
275 if (aSupportShape.IsSame(anExp.Current()))
281 // Construct a projection vertex.
282 const gp_Pnt &aPntProj = aDistShSh.PointOnShape2(i);
283 // TopoDS_Shape aProj = BRepBuilderAPI_MakeVertex(aPntProj).Shape();
284 DEBTRACE("projection: (" << aPntProj.X() << ", " << aPntProj.Y() << ", " << aPntProj.Z() << ")");
285 gp_XY aProjXY = gp_XY(aPntProj.X(), aPntProj.Y());
286 aProjXY.Subtract(thePoint);
287 double distance = aProjXY.Modulus();
288 DEBTRACE("distance " << distance);
292 HYDROData_ProfileUZ::PointsList aProfilePoints = aProfile->GetParametricPoints();
293 for ( int i = 1, aNbPoints = aProfilePoints.Size(); i <= aNbPoints; ++i )
295 const HYDROData_IPolyline::Point& aPolylinePoint = aProfilePoints.Value( i );
296 DEBTRACE(" profile point: " << aPolylinePoint.X() << " " << aPolylinePoint.Y());
297 if (aPolylinePoint.X() < distance)
300 pt1 = aPolylinePoint;
302 if (aPolylinePoint.X() >= distance)
304 pt2 = aPolylinePoint;
308 if ((i1 == 0) && (distance > 0))
310 DEBTRACE("Projection aborted : non centered profile");
314 // get the point with the same param on altitude profile
315 double zRef = 0; //aPntProj.Z();
316 TopoDS_Shape aShapeUZ = aGuideUZ->GetShape();
317 if (!aShapeUZ.IsNull())
319 TopExp_Explorer anExpUZ(aShapeUZ, TopAbs_EDGE);
320 TopoDS_Edge anEdgeUZ;
322 for (; anExpUZ.More(); anExpUZ.Next(), anI++)
326 anEdgeUZ = TopoDS::Edge(anExpUZ.Current());
330 if (!anEdgeUZ.IsNull())
332 double anEdgePars[2];
333 Handle(Geom_Curve) aCurve = BRep_Tool::Curve(anEdgeUZ, anEdgePars[0], anEdgePars[1]);
335 aCurve->D0(aParam, aVal);
336 DEBTRACE("altitude point: " << aVal.X() << " "<< aVal.Y() << " "<< aVal.Z());
343 // distance on section to get altitude
344 if (i1 == aProfilePoints.Size()) // distance >= profile width
346 aResAltitude = zRef + pt1.Y();
350 double z = pt1.Y() + (pt2.Y() - pt1.Y())*(distance -pt1.X())/(pt2.X()-pt1.X());
351 aResAltitude = zRef + z;
353 DEBTRACE("altitude: " << aResAltitude << " delta: " << aResAltitude - zRef);
358 if (!hasValidSolution)
360 DEBTRACE("Projection aborted : no projection");