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_HANDLE(HYDROData_ChannelAltitude, HYDROData_IAltitudeObject)
  57 IMPLEMENT_STANDARD_RTTIEXT(HYDROData_ChannelAltitude, HYDROData_IAltitudeObject)
  58
  59 HYDROData_ChannelAltitude::HYDROData_ChannelAltitude()
  60 : HYDROData_IAltitudeObject()
  61 {
  62 }
  63
  64 HYDROData_ChannelAltitude::~HYDROData_ChannelAltitude()
  65 {
  66 }
  67
  68 double HYDROData_ChannelAltitude::GetAltitudeForPoint( const gp_XY& thePoint ) const
  69 {
  70   DEBTRACE("HYDROData_ChannelAltitude::GetAltitudeForPoint");
  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   Handle(HYDROData_Polyline3D) aGuideLine = aChannel->GetGuideLine();
  82   if (aGuideLine.IsNull())
  83     {
  84       DEBTRACE("aGuideLine.IsNull()");
  85       return aResAltitude;
  86     }
  87   DEBTRACE("aGuideLine: " << aGuideLine->GetName().toStdString());
  88   Handle (HYDROData_Profile) aProfile = aChannel->GetProfile();
  89   if (aProfile.IsNull())
  90     {
  91       return aResAltitude;
  92     }
  93   DEBTRACE("aProfile: " << aProfile->GetName().toStdString());
  94
  95   // --- See GEOMImpl_ProjectionDriver.cxx
  96
  97   TopoDS_Shape aShape = aGuideLine->GetShape3D();
  98   double middleZ = -9999;
  99   aGuideLine->GetMiddleZ(middleZ); // use the middle Z value of the 3d line to help the projection.
 100   if (middleZ < -9000)
 101     {
 102       DEBTRACE("the middle Z value of the 3d line is incorrect");
 103     }
 104   gp_Pnt P1(thePoint.X(), thePoint.Y(), middleZ);
 105   TopoDS_Shape aPoint = BRepBuilderAPI_MakeVertex(P1).Shape();
 106
 107   if (aPoint.IsNull() || aShape.IsNull())
 108     {
 109       DEBTRACE("aPoint.IsNull() || aShape.IsNull()");
 110       return aResAltitude;
 111     }
 112
 113   if (aShape.ShapeType() != TopAbs_EDGE && aShape.ShapeType() != TopAbs_WIRE)
 114     {
 115       DEBTRACE("Projection aborted : the shape is neither an edge nor a wire");
 116       return aResAltitude;
 117     }
 118
 119   // Perform projection.
 120   BRepExtrema_DistShapeShape aDistShSh(aPoint, aShape, Extrema_ExtFlag_MIN);
 121
 122   if (aDistShSh.IsDone() == Standard_False)
 123     {
 124       DEBTRACE("Projection not done");
 125       return aResAltitude;
 126     }
 127
 128   Standard_Boolean hasValidSolution = Standard_False;
 129   Standard_Integer aNbSolutions = aDistShSh.NbSolution();
 130   Standard_Integer i;
 131   double aParam = 0.;
 132   Standard_Real aTolConf = BRep_Tool::Tolerance(TopoDS::Vertex(aPoint));
 133   Standard_Real aTolAng = 1.e-4;
 134
 135   for (i = 1; i <= aNbSolutions; i++)
 136     {
 137       Standard_Boolean isValid = Standard_False;
 138       BRepExtrema_SupportType aSupportType = aDistShSh.SupportTypeShape2(i);
 139       TopoDS_Shape aSupportShape = aDistShSh.SupportOnShape2(i);
 140
 141       if (aSupportType == BRepExtrema_IsOnEdge)
 142         {
 143           // Minimal distance inside edge is really a projection.
 144           isValid = Standard_True;
 145           aDistShSh.ParOnEdgeS2(i, aParam);
 146         }
 147       else if (aSupportType == BRepExtrema_IsVertex)
 148         {
 149           TopExp_Explorer anExp(aShape, TopAbs_EDGE);
 150
 151           if (aDistShSh.Value() <= aTolConf)
 152             {
 153               // The point lies on the shape. This means this point
 154               // is really a projection.
 155               for (; anExp.More() && !isValid; anExp.Next())
 156                 {
 157                   TopoDS_Edge aCurEdge = TopoDS::Edge(anExp.Current());
 158
 159                   if (aCurEdge.IsNull() == Standard_False)
 160                     {
 161                       TopoDS_Vertex aVtx[2];
 162
 163                       TopExp::Vertices(aCurEdge, aVtx[0], aVtx[1]);
 164
 165                       for (int j = 0; j < 2; j++)
 166                         {
 167                           if (aSupportShape.IsSame(aVtx[j]))
 168                             {
 169                               // The current edge is a projection edge.
 170                               isValid = Standard_True;
 171                               aSupportShape = aCurEdge;
 172                               aParam = BRep_Tool::Parameter(aVtx[j], aCurEdge);
 173                               break;
 174                             }
 175                         }
 176                     }
 177                 }
 178             }
 179           else
 180             {
 181               // Minimal distance to vertex is not always a real projection.
 182               gp_Pnt aPnt = BRep_Tool::Pnt(TopoDS::Vertex(aPoint));
 183               gp_Pnt aPrjPnt = BRep_Tool::Pnt(TopoDS::Vertex(aSupportShape));
 184               gp_Vec aDProjP(aPrjPnt, aPnt);
 185
 186               for (; anExp.More() && !isValid; anExp.Next())
 187                 {
 188                   TopoDS_Edge aCurEdge = TopoDS::Edge(anExp.Current());
 189
 190                   if (aCurEdge.IsNull() == Standard_False)
 191                     {
 192                       TopoDS_Vertex aVtx[2];
 193
 194                       TopExp::Vertices(aCurEdge, aVtx[0], aVtx[1]);
 195
 196                       for (int j = 0; j < 2; j++)
 197                         {
 198                           if (aSupportShape.IsSame(aVtx[j]))
 199                             {
 200                               // Check if the point is a projection to the current edge.
 201                               Standard_Real anEdgePars[2];
 202                               Handle(Geom_Curve) aCurve = BRep_Tool::Curve(aCurEdge, anEdgePars[0], anEdgePars[1]);
 203                               gp_Pnt aVal;
 204                               gp_Vec aD1;
 205
 206                               aParam = BRep_Tool::Parameter(aVtx[j], aCurEdge);
 207                               aCurve->D1(aParam, aVal, aD1);
 208
 209                               if (Abs(aD1.Dot(aDProjP)) <= aTolAng)
 210                                 {
 211                                   // The current edge is a projection edge.
 212                                   isValid = Standard_True;
 213                                   aSupportShape = aCurEdge;
 214                                   break;
 215                                 }
 216                             }
 217                         }
 218                     }
 219                 }
 220             }
 221         }
 222
 223       if (isValid)
 224         {
 225           if (hasValidSolution)
 226             {
 227               DEBTRACE("Projection aborted : multiple solutions");
 228               return aResAltitude;
 229             }
 230
 231           // Store the valid solution.
 232           hasValidSolution = Standard_True;
 233
 234           // Normalize parameter.
 235           TopoDS_Edge aSupportEdge = TopoDS::Edge(aSupportShape);
 236           Standard_Real aF, aL;
 237
 238           BRep_Tool::Range(aSupportEdge, aF, aL);
 239
 240           if (Abs(aL - aF) <= aTolConf)
 241             {
 242               DEBTRACE("Projection aborted : degenerated projection edge");
 243               return aResAltitude;
 244             }
 245
 246           // Construct a projection vertex.
 247           const gp_Pnt &aPntProj = aDistShSh.PointOnShape2(i);
 248           TopoDS_Shape aProj = BRepBuilderAPI_MakeVertex(aPntProj).Shape();
 249           DEBTRACE("projection: (" << aPntProj.X() << ", " << aPntProj.Y() << ", " << aPntProj.Z() << ")");
 250           gp_XY aProjXY = gp_XY(aPntProj.X(), aPntProj.Y());
 251           aProjXY.Subtract(thePoint);
 252           double distance = aProjXY.Modulus();
 253           DEBTRACE("distance " << distance);
 254           int i1 = 0;
 255           gp_XY pt1 = gp_XY();
 256           gp_XY pt2 = gp_XY();
 257           HYDROData_ProfileUZ::PointsList aProfilePoints = aProfile->GetParametricPoints();
 258           for ( int i = 1, aNbPoints = aProfilePoints.Size(); i <= aNbPoints; ++i )
 259             {
 260               const HYDROData_IPolyline::Point& aPolylinePoint = aProfilePoints.Value( i );
 261               DEBTRACE("  profile point: " << aPolylinePoint.X() << " " << aPolylinePoint.Y());
 262               if (aPolylinePoint.X() < distance)
 263                 {
 264                   i1 = i;
 265                   pt1 = aPolylinePoint;
 266                 }
 267               if (aPolylinePoint.X() >= distance)
 268                 {
 269                   pt2 = aPolylinePoint;
 270                   break;
 271                 }
 272             }
 273           if ((i1 == 0) && (distance > 0))
 274             {
 275               DEBTRACE("Projection aborted : non centered profile");
 276               return aResAltitude;
 277             }
 278           if (i1 == aProfilePoints.Size()) // distance >= profile width
 279             {
 280               aResAltitude = aPntProj.Z() + pt1.Y();
 281             }
 282           else
 283             {
 284               double z = pt1.Y() + (pt2.Y() - pt1.Y())*(distance -pt1.X())/(pt2.X()-pt1.X());
 285               aResAltitude = aPntProj.Z() + z;
 286             }
 287           DEBTRACE("altitude: " << aResAltitude << " delta: " << aResAltitude - aPntProj.Z());
 288           return aResAltitude;
 289         }
 290     }
 291
 292   if (!hasValidSolution)
 293     {
 294       DEBTRACE("Projection aborted : no projection");
 295       return aResAltitude;
 296     }
 297
 298   return aResAltitude;
 299 }
 300
 301
 302
 303