src/GEOMAlgo/GEOMAlgo_FinderShapeOnQuad.cxx

   1 // Copyright (C) 2007-2023  CEA, EDF, OPEN CASCADE
   2 //
   3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   6 // This library is free software; you can redistribute it and/or
   7 // modify it under the terms of the GNU Lesser General Public
   8 // License as published by the Free Software Foundation; either
   9 // version 2.1 of the License, or (at your option) any later version.
  10 //
  11 // This library is distributed in the hope that it will be useful,
  12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 // Lesser General Public License for more details.
  15 //
  16 // You should have received a copy of the GNU Lesser General Public
  17 // License along with this library; if not, write to the Free Software
  18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22
  23 // File      : GEOMAlgo_FinderShapeOnQuad.cxx
  24 // Created   : Mon Oct 17 17:31:45 2005
  25 // Author    : Edward AGAPOV (eap)
  26 //
  27 #include "GEOMAlgo_FinderShapeOnQuad.hxx"
  28 #include "GEOMAlgo_SurfaceTools.hxx"
  29
  30 #include <gp_Pnt.hxx>
  31 #include <gp_Vec.hxx>
  32 #include <Geom_Plane.hxx>
  33
  34
  35 GEOMAlgo_FinderShapeOnQuad::GEOMAlgo_FinderShapeOnQuad(const gp_Pnt & theTopLeftPoint,
  36                                                        const gp_Pnt & theTopRightPoint,
  37                                                        const gp_Pnt & theBottomLeftPoint,
  38                                                        const gp_Pnt & theBottomRightPoint)
  39 {
  40   myPoints.resize(6);
  41   myPoints[0] = theTopLeftPoint    ;
  42   myPoints[1] = theTopRightPoint   ;
  43   myPoints[2] = theBottomRightPoint;
  44   myPoints[3] = theBottomLeftPoint ;
  45   myPoints[4] = myPoints[0];
  46   myPoints[5] = myPoints[1];
  47
  48   // Find plane normal defined by corner points, it will be used to define a plane
  49   // for each quadrangle side.
  50   myQuadNormal.SetCoord (0,0,0);
  51   for ( int i = 1; i <= 4; ++i )
  52     myQuadNormal += gp_Vec( myPoints[i], myPoints[i+1] ) ^ gp_Vec( myPoints[i], myPoints[i-1] );
  53   //std::cout<<std::endl<<" X Vec : "<<myQuadNormal.X()<<" "<<myQuadNormal.Y()<<" "<<myQuadNormal.Z()<<" "<<endl;
  54
  55   if ( myQuadNormal.SquareMagnitude() <= DBL_MIN ) {
  56     myErrorStatus = 101;
  57     return;
  58   }
  59
  60   // detect concave quadrangle sides
  61   myConcaveQuad = false;
  62   myConcaveSide.resize (4, false);
  63   for ( int i = 1; i <= 4; ++i ) {
  64     gp_Vec localQN = gp_Vec( myPoints[i], myPoints[i+1] ) ^ gp_Vec( myPoints[i], myPoints[i-1] );
  65     if ( myQuadNormal * localQN < 0 )
  66       myConcaveSide[i-1] = myConcaveSide[i] = myConcaveQuad = true;
  67   }
  68
  69   // loop on quadrangle sides
  70   myPlanes.reserve( 4 );
  71   for ( int i = 0; i < 4; ++i )
  72   {
  73     // point1 -> point2 vector
  74     gp_Vec aSideVec( myPoints[ i ], myPoints[ i + 1 ]);
  75     //std::cout<<" Y Vec : "<<aSideVec.X()<<" "<<aSideVec.Y()<<" "<<aSideVec.Z()<<" "<<endl;
  76
  77     // plane normal
  78     gp_Vec aSideNorm = aSideVec ^ myQuadNormal;
  79     if ( aSideNorm.SquareMagnitude() <= DBL_MIN )
  80       continue;
  81     //std::cout<<" Z Vec : "<<aSideNorm.X()<<" "<<aSideNorm.Y()<<" "<<aSideNorm.Z()<<" "<<endl;
  82
  83     // make plane
  84     Handle(Geom_Plane) aPlane = new Geom_Plane( myPoints[ i ], aSideNorm );
  85     myPlanes.push_back( GeomAdaptor_Surface() );
  86     myPlanes.back().Load( aPlane );
  87   }
  88 }
  89
  90 //=======================================================================
  91 //function : CheckData
  92 //purpose  :
  93 //=======================================================================
  94
  95 void GEOMAlgo_FinderShapeOnQuad::CheckData()
  96 {
  97   if ( !myPlanes.empty() )
  98     mySurface = myPlanes[0].Surface();
  99   GEOMAlgo_FinderShapeOn1::CheckData();
 100 }
 101
 102 //=======================================================================
 103 //function : GetPointState
 104 //purpose  :
 105 //=======================================================================
 106
 107 TopAbs_State GEOMAlgo_FinderShapeOnQuad::GetPointState(const gp_Pnt& aP)
 108 {
 109   // Return IN if aP has TopAbs_IN with all sides.
 110   // In the case of concave quadrangle, return IN if
 111   // aP is OUT of only one concave side
 112   double nbIn = 0.;
 113   for ( size_t i = 0; i < myPlanes.size(); ++i )
 114   {
 115     TopAbs_State aSt;
 116     GEOMAlgo_SurfaceTools::GetState(aP, myPlanes[i], myTolerance, aSt);
 117     if ( aSt == TopAbs_IN )
 118     {
 119       nbIn += myConcaveSide[i] ? 0.5 : 1.0;
 120     }
 121     else if ( aSt == TopAbs_ON )
 122     {
 123       // check that aP is between quadrangle corners
 124       Handle(Geom_Plane) aSidePlane = Handle(Geom_Plane)::DownCast( myPlanes[i].Surface() );
 125       gp_Vec aSideNorm = aSidePlane->Axis().Direction();
 126       gp_Vec aSideVec = myQuadNormal ^ aSideNorm;
 127       gp_Vec c1p ( myPoints[i], aP );
 128       gp_Vec pc2 ( aP, myPoints[i+1] );
 129       if ( aSideVec * c1p >= 0. && aSideVec * pc2 >= 0. )
 130         return TopAbs_ON;
 131       // consider to be IN (???????????)
 132       //nbIn += myConcaveSide[i] ? 0.5 : 1.0;
 133     }
 134   }
 135   Standard_Real inThreshold = myPlanes.size(); // usually 4.0
 136   if ( myConcaveQuad )
 137     inThreshold = 2.5; // 1.0 + 1.0 + 0.5
 138
 139   if ( nbIn >= inThreshold )
 140     return TopAbs_IN;
 141
 142   return TopAbs_OUT;
 143 }
 144