Version: 6.3.1

src/SMESH/SMESH_Block.cxx

Go to the documentation of this file.
00001 // Copyright (C) 2007-2011  CEA/DEN, EDF R&D, OPEN CASCADE
00002 //
00003 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
00004 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
00005 //
00006 // This library is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 2.1 of the License.
00010 //
00011 // This library is distributed in the hope that it will be useful,
00012 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00013 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00014 // Lesser General Public License for more details.
00015 //
00016 // You should have received a copy of the GNU Lesser General Public
00017 // License along with this library; if not, write to the Free Software
00018 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00019 //
00020 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00021 //
00022 
00023 // File      : SMESH_Pattern.hxx
00024 // Created   : Mon Aug  2 10:30:00 2004
00025 // Author    : Edward AGAPOV (eap)
00026 //
00027 #include "SMESH_Block.hxx"
00028 
00029 #include <BRepAdaptor_Curve.hxx>
00030 #include <BRepAdaptor_Curve2d.hxx>
00031 #include <BRepAdaptor_Surface.hxx>
00032 #include <BRepTools.hxx>
00033 #include <BRepTools_WireExplorer.hxx>
00034 #include <BRep_Builder.hxx>
00035 #include <BRep_Tool.hxx>
00036 #include <Bnd_Box.hxx>
00037 #include <Extrema_ExtPC.hxx>
00038 #include <Extrema_POnCurv.hxx>
00039 #include <Geom2d_Curve.hxx>
00040 #include <ShapeAnalysis.hxx>
00041 #include <TopExp_Explorer.hxx>
00042 #include <TopTools_IndexedDataMapOfShapeListOfShape.hxx>
00043 #include <TopTools_ListIteratorOfListOfShape.hxx>
00044 #include <TopTools_ListOfShape.hxx>
00045 #include <TopoDS.hxx>
00046 #include <TopoDS_Compound.hxx>
00047 #include <TopoDS_Face.hxx>
00048 #include <TopoDS_Iterator.hxx>
00049 #include <TopoDS_Wire.hxx>
00050 #include <gp_Trsf.hxx>
00051 #include <gp_Vec.hxx>
00052 #include <math_FunctionSetRoot.hxx>
00053 #include <math_Matrix.hxx>
00054 #include <math_Vector.hxx>
00055 
00056 #include "SMDS_MeshNode.hxx"
00057 #include "SMDS_MeshVolume.hxx"
00058 #include "SMDS_VolumeTool.hxx"
00059 #include "utilities.h"
00060 
00061 #include <list>
00062 
00063 using namespace std;
00064 
00065 //#define DEBUG_PARAM_COMPUTE
00066 
00067 //================================================================================
00074 //================================================================================
00075 
00076 void SMESH_Block::TEdge::Set( const int edgeID, Adaptor3d_Curve* curve, const bool isForward )
00077 {
00078   myCoordInd = SMESH_Block::GetCoordIndOnEdge( edgeID );
00079   if ( myC3d ) delete myC3d;
00080   myC3d = curve;
00081   myFirst = curve->FirstParameter();
00082   myLast = curve->LastParameter();
00083   if ( !isForward )
00084     std::swap( myFirst, myLast );
00085 }
00086 
00087 //================================================================================
00094 //================================================================================
00095 
00096 void SMESH_Block::TEdge::Set( const int edgeID, const gp_XYZ& node1, const gp_XYZ& node2 )
00097 {
00098   myCoordInd = SMESH_Block::GetCoordIndOnEdge( edgeID );
00099   myNodes[ 0 ] = node1;
00100   myNodes[ 1 ] = node2;
00101 
00102   if ( myC3d ) delete myC3d;
00103   myC3d = 0;
00104 }
00105 
00106 //=======================================================================
00107 //function : SMESH_Block::TEdge::GetU
00108 //purpose  : 
00109 //=======================================================================
00110 
00111 double SMESH_Block::TEdge::GetU( const gp_XYZ& theParams ) const
00112 {
00113   double u = theParams.Coord( myCoordInd );
00114   if ( !myC3d ) // if mesh block
00115     return u;
00116   return ( 1 - u ) * myFirst + u * myLast;
00117 }
00118 
00119 //=======================================================================
00120 //function : SMESH_Block::TEdge::Point
00121 //purpose  : 
00122 //=======================================================================
00123 
00124 gp_XYZ SMESH_Block::TEdge::Point( const gp_XYZ& theParams ) const
00125 {
00126   double u = GetU( theParams );
00127   if ( myC3d ) return myC3d->Value( u ).XYZ();
00128   // mesh block
00129   return myNodes[0] * ( 1 - u ) + myNodes[1] * u;
00130 }
00131 
00132 //================================================================================
00136 //================================================================================
00137 
00138 SMESH_Block::TEdge::~TEdge()
00139 {
00140   if ( myC3d ) delete myC3d;
00141 }
00142 
00143 //================================================================================
00151 //================================================================================
00152 
00153 void SMESH_Block::TFace::Set( const int          faceID,
00154                               Adaptor3d_Surface* S,
00155                               Adaptor2d_Curve2d* c2D[4],
00156                               const bool         isForward[4] )
00157 {
00158   if ( myS ) delete myS;
00159   myS = S;
00160   // pcurves
00161   vector< int > edgeIdVec;
00162   GetFaceEdgesIDs( faceID, edgeIdVec );
00163   for ( int iE = 0; iE < edgeIdVec.size(); iE++ ) // loop on 4 edges
00164   {
00165     myCoordInd[ iE ] = GetCoordIndOnEdge( edgeIdVec[ iE ] );
00166     if ( myC2d[ iE ]) delete myC2d[ iE ];
00167     myC2d[ iE ] = c2D[ iE ];
00168     myFirst[ iE ] = myC2d[ iE ]->FirstParameter();
00169     myLast [ iE ] = myC2d[ iE ]->LastParameter();
00170     if ( !isForward[ iE ])
00171       std::swap( myFirst[ iE ], myLast[ iE ] );
00172   }
00173   // 2d corners
00174   myCorner[ 0 ] = myC2d[ 0 ]->Value( myFirst[0] ).XY();
00175   myCorner[ 1 ] = myC2d[ 0 ]->Value( myLast[0] ).XY();
00176   myCorner[ 2 ] = myC2d[ 1 ]->Value( myLast[1] ).XY();
00177   myCorner[ 3 ] = myC2d[ 1 ]->Value( myFirst[1] ).XY();
00178 }
00179 
00180 //================================================================================
00187 //================================================================================
00188 
00189 void SMESH_Block::TFace::Set( const int faceID, const TEdge& edgeU0, const TEdge& edgeU1 )
00190 {
00191   vector< int > edgeIdVec;
00192   GetFaceEdgesIDs( faceID, edgeIdVec );
00193   myNodes[ 0 ] = edgeU0.NodeXYZ( 1 );
00194   myNodes[ 1 ] = edgeU0.NodeXYZ( 0 );
00195   myNodes[ 2 ] = edgeU1.NodeXYZ( 0 );
00196   myNodes[ 3 ] = edgeU1.NodeXYZ( 1 );
00197   myCoordInd[ 0 ] = GetCoordIndOnEdge( edgeIdVec[ 0 ] );
00198   myCoordInd[ 1 ] = GetCoordIndOnEdge( edgeIdVec[ 1 ] );
00199   myCoordInd[ 2 ] = GetCoordIndOnEdge( edgeIdVec[ 2 ] );
00200   myCoordInd[ 3 ] = GetCoordIndOnEdge( edgeIdVec[ 3 ] );
00201   if ( myS ) delete myS;
00202   myS = 0;
00203 }
00204 
00205 //================================================================================
00209 //================================================================================
00210 
00211 SMESH_Block::TFace::~TFace()
00212 {
00213   if ( myS ) delete myS;
00214   for ( int i = 0 ; i < 4; ++i )
00215     if ( myC2d[ i ]) delete myC2d[ i ];
00216 }
00217 
00218 //=======================================================================
00219 //function : SMESH_Block::TFace::GetCoefs
00220 //purpose  : return coefficients for addition of [0-3]-th edge and vertex
00221 //=======================================================================
00222 
00223 void SMESH_Block::TFace::GetCoefs(int           iE,
00224                                   const gp_XYZ& theParams,
00225                                   double&       Ecoef,
00226                                   double&       Vcoef ) const
00227 {
00228   double dU = theParams.Coord( GetUInd() );
00229   double dV = theParams.Coord( GetVInd() );
00230   switch ( iE ) {
00231   case 0:
00232     Ecoef = ( 1 - dV ); // u0
00233     Vcoef = ( 1 - dU ) * ( 1 - dV ); break; // 00
00234   case 1:
00235     Ecoef = dV; // u1
00236     Vcoef = dU * ( 1 - dV ); break; // 10
00237   case 2:
00238     Ecoef = ( 1 - dU ); // 0v
00239     Vcoef = dU * dV  ; break; // 11
00240   case 3:
00241     Ecoef = dU  ; // 1v
00242     Vcoef = ( 1 - dU ) * dV  ; break; // 01
00243   default: ASSERT(0);
00244   }
00245 }
00246 
00247 //=======================================================================
00248 //function : SMESH_Block::TFace::GetUV
00249 //purpose  : 
00250 //=======================================================================
00251 
00252 gp_XY SMESH_Block::TFace::GetUV( const gp_XYZ& theParams ) const
00253 {
00254   gp_XY uv(0.,0.);
00255   for ( int iE = 0; iE < 4; iE++ ) // loop on 4 edges
00256   {
00257     double Ecoef = 0, Vcoef = 0;
00258     GetCoefs( iE, theParams, Ecoef, Vcoef );
00259     // edge addition
00260     double u = theParams.Coord( myCoordInd[ iE ] );
00261     u = ( 1 - u ) * myFirst[ iE ] + u * myLast[ iE ];
00262     uv += Ecoef * myC2d[ iE ]->Value( u ).XY();
00263     // corner addition
00264     uv -= Vcoef * myCorner[ iE ];
00265   }
00266   return uv;
00267 }
00268 
00269 //=======================================================================
00270 //function : SMESH_Block::TFace::Point
00271 //purpose  : 
00272 //=======================================================================
00273 
00274 gp_XYZ SMESH_Block::TFace::Point( const gp_XYZ& theParams ) const
00275 {
00276   gp_XYZ p(0.,0.,0.);
00277   if ( !myS ) // if mesh block
00278   {
00279     for ( int iE = 0; iE < 4; iE++ ) // loop on 4 edges
00280     {
00281       double Ecoef = 0, Vcoef = 0;
00282       GetCoefs( iE, theParams, Ecoef, Vcoef );
00283       // edge addition
00284       double u = theParams.Coord( myCoordInd[ iE ] );
00285       int i1 = 0, i2 = 1;
00286       switch ( iE ) {
00287       case 1: i1 = 3; i2 = 2; break;
00288       case 2: i1 = 1; i2 = 2; break;
00289       case 3: i1 = 0; i2 = 3; break;
00290       }
00291       p += Ecoef * ( myNodes[ i1 ] * ( 1 - u ) + myNodes[ i2 ] * u );
00292       // corner addition
00293       p -= Vcoef * myNodes[ iE ];
00294     }
00295     
00296   }
00297   else // shape block
00298   {
00299     gp_XY uv = GetUV( theParams );
00300     p = myS->Value( uv.X(), uv.Y() ).XYZ();
00301   }
00302   return p;
00303 }
00304 
00305 //=======================================================================
00306 //function : GetShapeCoef
00307 //purpose  : 
00308 //=======================================================================
00309 
00310 double* SMESH_Block::GetShapeCoef (const int theShapeID)
00311 {
00312   static double shapeCoef[][3] = {
00313     //    V000,        V100,        V010,         V110
00314     { -1,-1,-1 }, {  1,-1,-1 }, { -1, 1,-1 }, {  1, 1,-1 },
00315     //    V001,        V101,        V011,         V111,
00316     { -1,-1, 1 }, {  1,-1, 1 }, { -1, 1, 1 }, {  1, 1, 1 },
00317     //    Ex00,        Ex10,        Ex01,         Ex11,
00318     {  0,-1,-1 }, {  0, 1,-1 }, {  0,-1, 1 }, {  0, 1, 1 },
00319     //    E0y0,        E1y0,        E0y1,         E1y1,
00320     { -1, 0,-1 }, {  1, 0,-1 }, { -1, 0, 1 }, {  1, 0, 1 },
00321     //    E00z,        E10z,        E01z,         E11z,
00322     { -1,-1, 0 }, {  1,-1, 0 }, { -1, 1, 0 }, {  1, 1, 0 },
00323     //    Fxy0,        Fxy1,        Fx0z,         Fx1z,         F0yz,           F1yz,
00324     {  0, 0,-1 }, {  0, 0, 1 }, {  0,-1, 0 }, {  0, 1, 0 }, { -1, 0, 0 }, {  1, 0, 0 },
00325     // ID_Shell
00326     {  0, 0, 0 }
00327   };
00328   if ( theShapeID < ID_V000 || theShapeID > ID_F1yz )
00329     return shapeCoef[ ID_Shell - 1 ];
00330 
00331   return shapeCoef[ theShapeID - 1 ];
00332 }
00333 
00334 //=======================================================================
00335 //function : ShellPoint
00336 //purpose  : return coordinates of a point in shell
00337 //=======================================================================
00338 
00339 bool SMESH_Block::ShellPoint( const gp_XYZ& theParams, gp_XYZ& thePoint ) const
00340 {
00341   thePoint.SetCoord( 0., 0., 0. );
00342   for ( int shapeID = ID_V000; shapeID < ID_Shell; shapeID++ )
00343   {
00344     // coef
00345     double* coefs = GetShapeCoef( shapeID );
00346     double k = 1;
00347     for ( int iCoef = 0; iCoef < 3; iCoef++ ) {
00348       if ( coefs[ iCoef ] != 0 ) {
00349         if ( coefs[ iCoef ] < 0 )
00350           k *= ( 1. - theParams.Coord( iCoef + 1 ));
00351         else
00352           k *= theParams.Coord( iCoef + 1 );
00353       }
00354     }
00355     // add point on a shape
00356     if ( fabs( k ) > DBL_MIN )
00357     {
00358       gp_XYZ Ps;
00359       if ( shapeID < ID_Ex00 ) // vertex
00360         VertexPoint( shapeID, Ps );
00361       else if ( shapeID < ID_Fxy0 ) { // edge
00362         EdgePoint( shapeID, theParams, Ps );
00363         k = -k;
00364       } else // face
00365         FacePoint( shapeID, theParams, Ps );
00366 
00367       thePoint += k * Ps;
00368     }
00369   }
00370   return true;
00371 }
00372 
00373 //=======================================================================
00374 //function : ShellPoint
00375 //purpose  : computes coordinates of a point in shell by points on sub-shapes;
00376 //           thePointOnShape[ subShapeID ] must be a point on a subShape
00377 //=======================================================================
00378 
00379 bool SMESH_Block::ShellPoint(const gp_XYZ&         theParams,
00380                              const vector<gp_XYZ>& thePointOnShape,
00381                              gp_XYZ&               thePoint )
00382 {
00383   if ( thePointOnShape.size() < ID_F1yz )
00384     return false;
00385 
00386   const double x = theParams.X(), y = theParams.Y(), z = theParams.Z();
00387   const double x1 = 1. - x,       y1 = 1. - y,       z1 = 1. - z;
00388   const vector<gp_XYZ>& p = thePointOnShape;
00389 
00390   thePoint = 
00391     x1 * p[ID_F0yz] + x * p[ID_F1yz] +
00392     y1 * p[ID_Fx0z] + y * p[ID_Fx1z] +
00393     z1 * p[ID_Fxy0] + z * p[ID_Fxy1] +
00394     x1 * (y1 * (z1 * p[ID_V000] + z * p[ID_V001])  +
00395           y  * (z1 * p[ID_V010] + z * p[ID_V011])) +
00396     x  * (y1 * (z1 * p[ID_V100] + z * p[ID_V101])  +
00397           y  * (z1 * p[ID_V110] + z * p[ID_V111]));
00398   thePoint -=
00399     x1 * (y1 * p[ID_E00z] + y * p[ID_E01z]) + 
00400     x  * (y1 * p[ID_E10z] + y * p[ID_E11z]) + 
00401     y1 * (z1 * p[ID_Ex00] + z * p[ID_Ex01]) + 
00402     y  * (z1 * p[ID_Ex10] + z * p[ID_Ex11]) + 
00403     z1 * (x1 * p[ID_E0y0] + x * p[ID_E1y0]) + 
00404     z  * (x1 * p[ID_E0y1] + x * p[ID_E1y1]);
00405 
00406   return true;
00407 }
00408 
00409 //=======================================================================
00410 //function : Constructor
00411 //purpose  : 
00412 //=======================================================================
00413 
00414 SMESH_Block::SMESH_Block():
00415   myNbIterations(0),
00416   mySumDist(0.),
00417   myTolerance(-1.) // to be re-initialized
00418 {
00419 }
00420 
00421 
00422 //=======================================================================
00423 //function : NbVariables
00424 //purpose  : 
00425 //=======================================================================
00426 
00427 Standard_Integer SMESH_Block::NbVariables() const
00428 {
00429   return 3;
00430 }
00431 
00432 //=======================================================================
00433 //function : NbEquations
00434 //purpose  : 
00435 //=======================================================================
00436 
00437 Standard_Integer SMESH_Block::NbEquations() const
00438 {
00439   return 1;
00440 }
00441 
00442 //=======================================================================
00443 //function : Value
00444 //purpose  : 
00445 //=======================================================================
00446 
00447 Standard_Boolean SMESH_Block::Value(const math_Vector& theXYZ, math_Vector& theFxyz) 
00448 {
00449   gp_XYZ P, params( theXYZ(1), theXYZ(2), theXYZ(3) );
00450   if ( params.IsEqual( myParam, DBL_MIN )) { // same param
00451     theFxyz( 1 ) = funcValue( myValues[ SQUARE_DIST ]);
00452   }
00453   else {
00454     ShellPoint( params, P );
00455     gp_Vec dP( P - myPoint );
00456     theFxyz(1) = funcValue( dP.SquareMagnitude() );
00457   }
00458   return true;
00459 }
00460 
00461 //=======================================================================
00462 //function : Derivatives
00463 //purpose  : 
00464 //=======================================================================
00465 
00466 Standard_Boolean SMESH_Block::Derivatives(const math_Vector& XYZ,math_Matrix& Df) 
00467 {
00468   math_Vector F(1,3);
00469   return Values(XYZ,F,Df);
00470 }
00471 
00472 //=======================================================================
00473 //function : GetStateNumber
00474 //purpose  : 
00475 //=======================================================================
00476 
00477 Standard_Integer SMESH_Block::GetStateNumber ()
00478 {
00479   return 0; //myValues[0] < 1e-1;
00480 }
00481 
00482 //=======================================================================
00483 //function : Values
00484 //purpose  : 
00485 //=======================================================================
00486 
00487 Standard_Boolean SMESH_Block::Values(const math_Vector& theXYZ,
00488                                      math_Vector&       theFxyz,
00489                                      math_Matrix&       theDf) 
00490 {
00491   gp_XYZ P, params( theXYZ(1), theXYZ(2), theXYZ(3) );
00492   if ( params.IsEqual( myParam, DBL_MIN )) { // same param
00493     theFxyz( 1 )      = funcValue( myValues[ SQUARE_DIST ] );
00494     theDf( 1, DRV_1 ) = myValues[ DRV_1 ];
00495     theDf( 1, DRV_2 ) = myValues[ DRV_2 ];
00496     theDf( 1, DRV_3 ) = myValues[ DRV_3 ];
00497     return true;
00498   }
00499 #ifdef DEBUG_PARAM_COMPUTE
00500   MESSAGE ( "PARAM GUESS: " << params.X() << " "<< params.Y() << " "<< params.X() );
00501   myNbIterations++; // how many times call ShellPoint()
00502 #endif
00503   ShellPoint( params, P );
00504 
00505   gp_Vec dP( myPoint, P );
00506   double sqDist = dP.SquareMagnitude();
00507   theFxyz(1) = funcValue( sqDist );
00508 
00509   if ( sqDist < myTolerance * myTolerance ) { // a solution found
00510     myParam = params;
00511     myValues[ SQUARE_DIST ] = sqDist;
00512     theFxyz(1)  = theDf( 1,1 ) = theDf( 1,2 ) = theDf( 1,3 ) = 0;
00513     return true;
00514   }
00515 
00516   if ( sqDist < myValues[ SQUARE_DIST ] ) // a better guess
00517   {
00518     // 3 partial derivatives
00519     gp_Vec drv[ 3 ]; // where we move with a small step in each direction
00520     for ( int iP = 1; iP <= 3; iP++ ) {
00521       if ( iP == myFaceIndex ) {
00522         drv[ iP - 1 ] = gp_Vec(0,0,0);
00523         continue;
00524       }
00525       gp_XYZ Pi;
00526       bool onEdge = ( theXYZ( iP ) + 0.001 > 1. );
00527       if ( onEdge )
00528         params.SetCoord( iP, theXYZ( iP ) - 0.001 );
00529       else
00530         params.SetCoord( iP, theXYZ( iP ) + 0.001 );
00531       ShellPoint( params, Pi );
00532       params.SetCoord( iP, theXYZ( iP ) ); // restore params
00533       gp_Vec dPi ( P, Pi );
00534       if ( onEdge ) dPi *= -1.;
00535       double mag = dPi.Magnitude();
00536       if ( mag > DBL_MIN )
00537         dPi /= mag;
00538       drv[ iP - 1 ] = dPi;
00539     }
00540     for ( int iP = 0; iP < 3; iP++ ) {
00541 #if 1
00542       theDf( 1, iP + 1 ) = dP * drv[iP];
00543 #else
00544       // Distance from P to plane passing through myPoint and defined
00545       // by the 2 other derivative directions:
00546       // like IntAna_IntConicQuad::Perform (const gp_Lin& L, const gp_Pln& P)
00547       // where L is (P -> myPoint), P is defined by the 2 other derivative direction
00548       int iPrev = ( iP ? iP - 1 : 2 );
00549       int iNext = ( iP == 2 ? 0 : iP + 1 );
00550       gp_Vec plnNorm = drv[ iPrev ].Crossed( drv [ iNext ] );
00551       double Direc = plnNorm * drv[ iP ];
00552       if ( Abs(Direc) <= DBL_MIN )
00553         theDf( 1, iP + 1 ) = dP * drv[ iP ];
00554       else {
00555         double Dis = plnNorm * P - plnNorm * myPoint;
00556         theDf( 1, iP + 1 ) = Dis/Direc;
00557       }
00558 #endif
00559     }
00560 #ifdef DEBUG_PARAM_COMPUTE
00561     MESSAGE ( "F = " << theFxyz(1) << " DRV: " << theDf(1,1) << " " << theDf(1,2) << " " << theDf(1,3) );
00562     myNbIterations +=3; // how many times call ShellPoint()
00563 #endif
00564 
00565     // store better values
00566     myParam              = params;
00567     myValues[SQUARE_DIST]= sqDist;
00568     myValues[DRV_1]      = theDf(1,DRV_1);
00569     myValues[DRV_2]      = theDf(1,DRV_2);
00570     myValues[DRV_3]      = theDf(1,DRV_3);
00571   }
00572 
00573   return true;
00574 }
00575 
00576 //============================================================================
00577 //function : computeParameters
00578 //purpose  : compute point parameters in the block using math_FunctionSetRoot
00579 //============================================================================
00580 
00581 bool SMESH_Block::computeParameters(const gp_Pnt& thePoint,
00582                                     gp_XYZ&       theParams,
00583                                     const gp_XYZ& theParamsHint)
00584 {
00585   myPoint = thePoint.XYZ();
00586 
00587   myParam.SetCoord( -1,-1,-1 );
00588   myValues[ SQUARE_DIST ] = 1e100;
00589 
00590   math_Vector low  ( 1, 3, 0.0 );
00591   math_Vector up   ( 1, 3, 1.0 );
00592   math_Vector tol  ( 1, 3, 1e-4 );
00593   math_Vector start( 1, 3, 0.0 );
00594   start( 1 ) = theParamsHint.X();
00595   start( 2 ) = theParamsHint.Y();
00596   start( 3 ) = theParamsHint.Z();
00597 
00598   math_FunctionSetRoot paramSearch( *this, tol );
00599 
00600   mySquareFunc = 0; // large approaching steps
00601   //if ( hasHint ) mySquareFunc = 1; // small approaching steps
00602 
00603   double loopTol = 10 * myTolerance;
00604   int nbLoops = 0;
00605   while ( distance() > loopTol && nbLoops <= 3 )
00606   {
00607     paramSearch.Perform ( *static_cast<math_FunctionSetWithDerivatives*>(this),
00608                           start, low, up );
00609     start( 1 ) = myParam.X();
00610     start( 2 ) = myParam.Y();
00611     start( 3 ) = myParam.Z();
00612     mySquareFunc = !mySquareFunc;
00613     nbLoops++;
00614   }
00615 #ifdef DEBUG_PARAM_COMPUTE
00616   mySumDist += distance();
00617   MESSAGE ( " ------ SOLUTION: ( "<< myParam.X() <<" "<< myParam.Y() <<" "<< myParam.Z() <<" )"<<endl
00618          << " ------ DIST : " << distance() << "\t Tol=" << myTolerance << "\t Nb LOOPS=" << nbLoops << endl
00619          << " ------ NB IT: " << myNbIterations << ",  SUM DIST: " << mySumDist );
00620 #endif
00621 
00622   theParams = myParam;
00623 
00624   if ( myFaceIndex > 0 )
00625     theParams.SetCoord( myFaceIndex, myFaceParam );
00626 
00627   return true;
00628 }
00629 
00630 //=======================================================================
00631 //function : ComputeParameters
00632 //purpose  : compute point parameters in the block
00633 //=======================================================================
00634 
00635 bool SMESH_Block::ComputeParameters(const gp_Pnt& thePoint,
00636                                     gp_XYZ&       theParams,
00637                                     const int     theShapeID,
00638                                     const gp_XYZ& theParamsHint)
00639 {
00640   if ( VertexParameters( theShapeID, theParams ))
00641     return true;
00642 
00643   if ( IsEdgeID( theShapeID )) {
00644     TEdge& e = myEdge[ theShapeID - ID_FirstE ];
00645     Adaptor3d_Curve* curve = e.GetCurve();
00646     Extrema_ExtPC anExtPC( thePoint, *curve, curve->FirstParameter(), curve->LastParameter() );
00647     int i, nb = anExtPC.IsDone() ? anExtPC.NbExt() : 0;
00648     for ( i = 1; i <= nb; i++ ) {
00649       if ( anExtPC.IsMin( i ))
00650         return EdgeParameters( theShapeID, anExtPC.Point( i ).Parameter(), theParams );
00651     }
00652     return false;
00653   }
00654 
00655   const bool isOnFace = IsFaceID( theShapeID );
00656   double * coef = GetShapeCoef( theShapeID );
00657 
00658   // Find the first guess paremeters
00659 
00660   gp_XYZ start(0, 0, 0);
00661 
00662   bool hasHint = ( 0 <= theParamsHint.X() && theParamsHint.X() <= 1 &&
00663                    0 <= theParamsHint.Y() && theParamsHint.Y() <= 1 &&
00664                    0 <= theParamsHint.Y() && theParamsHint.Y() <= 1 );
00665   if ( !hasHint && !myGridComputed )
00666   {
00667     // define the first guess by thePoint projection on lines
00668     // connecting vertices
00669     bool needGrid = false;
00670     gp_XYZ par000( 0, 0, 0 ), par111( 1, 1, 1 );
00671     double zero = DBL_MIN * DBL_MIN;
00672     for ( int iEdge = 0, iParam = 1; iParam <= 3 && !needGrid; iParam++ )
00673     {
00674       if ( isOnFace && coef[ iParam - 1 ] != 0 ) {
00675         iEdge += 4;
00676         continue;
00677       }
00678       double sumParam = 0;
00679       for ( int iE = 0; iE < 4; iE++, iEdge++ ) { // loop on 4 parallel edges
00680         gp_Pnt p0 = myEdge[ iEdge ].Point( par000 );
00681         gp_Pnt p1 = myEdge[ iEdge ].Point( par111 );
00682         gp_Vec v01( p0, p1 ), v0P( p0, thePoint );
00683         double len2 = v01.SquareMagnitude();
00684         double par = 0;
00685         if ( len2 > zero ) {
00686           par = v0P.Dot( v01 ) / len2;
00687           if ( par < 0 || par > 1 ) { // projection falls out of line ends => needGrid
00688             needGrid = true;
00689             break;
00690           }
00691         }
00692         sumParam += par;
00693       }
00694       start.SetCoord( iParam, sumParam / 4.);
00695     }
00696     if ( needGrid ) {
00697       // compute nodes of 3 x 3 x 3 grid
00698       int iNode = 0;
00699       Bnd_Box box;
00700       for ( double x = 0.25; x < 0.9; x += 0.25 )
00701         for ( double y = 0.25; y < 0.9; y += 0.25 )
00702           for ( double z = 0.25; z < 0.9; z += 0.25 ) {
00703             TxyzPair & prmPtn = my3x3x3GridNodes[ iNode++ ];
00704             prmPtn.first.SetCoord( x, y, z );
00705             ShellPoint( prmPtn.first, prmPtn.second );
00706             box.Add( gp_Pnt( prmPtn.second ));
00707           }
00708       myGridComputed = true;
00709       myTolerance = sqrt( box.SquareExtent() ) * 1e-5;
00710     }
00711   }
00712 
00713   if ( hasHint )
00714   {
00715     start = theParamsHint;
00716   }
00717   else if ( myGridComputed )
00718   {
00719     double minDist = DBL_MAX;
00720     gp_XYZ* bestParam = 0;
00721     for ( int iNode = 0; iNode < 27; iNode++ ) {
00722       TxyzPair & prmPtn = my3x3x3GridNodes[ iNode ];
00723       double dist = ( thePoint.XYZ() - prmPtn.second ).SquareModulus();
00724       if ( dist < minDist ) {
00725         minDist = dist;
00726         bestParam = & prmPtn.first;
00727       }
00728     }
00729     start = *bestParam;
00730   }
00731 
00732   myFaceIndex = -1;
00733   myFaceParam = 0.;
00734   if ( isOnFace ) {
00735     // put a point on the face
00736     for ( int iCoord = 0; iCoord < 3; iCoord++ )
00737       if ( coef[ iCoord ] ) {
00738         myFaceIndex = iCoord + 1;
00739         myFaceParam = ( coef[ iCoord ] < 0.5 ) ? 0.0 : 1.0;
00740         start.SetCoord( myFaceIndex, myFaceParam );
00741       }
00742   }
00743 
00744 #ifdef DEBUG_PARAM_COMPUTE
00745   MESSAGE ( " #### POINT " <<thePoint.X()<<" "<<thePoint.Y()<<" "<<thePoint.Z()<<" ####" );
00746 #endif
00747 
00748   if ( myTolerance < 0 ) myTolerance = 1e-6;
00749 
00750   const double parDelta = 1e-4;
00751   const double sqTolerance = myTolerance * myTolerance;
00752 
00753   gp_XYZ solution = start, params = start;
00754   double sqDistance = 1e100; 
00755   int nbLoops = 0, nbGetWorst = 0;
00756 
00757   while ( nbLoops <= 100 )
00758   {
00759     gp_XYZ P, Pi;
00760     ShellPoint( params, P );
00761 
00762     gp_Vec dP( thePoint, P );
00763     double sqDist = dP.SquareMagnitude();
00764 
00765     if ( sqDist > sqDistance ) { // solution get worse
00766       if ( ++nbGetWorst > 2 )
00767         return computeParameters( thePoint, theParams, solution );
00768     }
00769 #ifdef DEBUG_PARAM_COMPUTE
00770     MESSAGE ( "PARAMS: ( " << params.X() <<" "<< params.Y() <<" "<< params.Z() <<" )" );
00771     MESSAGE ( "DIST: " << sqrt( sqDist ) );
00772 #endif
00773 
00774     if ( sqDist < sqDistance ) { // get better
00775       sqDistance = sqDist;
00776       solution   = params;
00777       nbGetWorst = 0;
00778       if ( sqDistance < sqTolerance ) // a solution found
00779         break;
00780     }
00781 
00782         // look for a next better solution
00783     for ( int iP = 1; iP <= 3; iP++ ) {
00784       if ( iP == myFaceIndex )
00785         continue;
00786       // see where we move with a small (=parDelta) step in this direction
00787       gp_XYZ nearParams = params;
00788       bool onEdge = ( params.Coord( iP ) + parDelta > 1. );
00789       if ( onEdge )
00790         nearParams.SetCoord( iP, params.Coord( iP ) - parDelta );
00791       else
00792         nearParams.SetCoord( iP, params.Coord( iP ) + parDelta );
00793       ShellPoint( nearParams, Pi );
00794       gp_Vec dPi ( P, Pi );
00795       if ( onEdge ) dPi *= -1.;
00796       // modify a parameter
00797       double mag = dPi.Magnitude();
00798       if ( mag < DBL_MIN )
00799         continue;
00800       gp_Vec dir = dPi / mag; // dir we move modifying the parameter
00801       double dist = dir * dP; // where we should get to
00802       double dPar = dist / mag * parDelta; // predict parameter change
00803       double curPar = params.Coord( iP );
00804       double par = curPar - dPar; // new parameter value
00805       while ( par > 1 || par < 0 ) {
00806         dPar /= 2.;
00807         par = curPar - dPar;
00808       }
00809       params.SetCoord( iP, par );
00810     }
00811 
00812     nbLoops++;
00813   }
00814 #ifdef DEBUG_PARAM_COMPUTE
00815   myNbIterations += nbLoops*4; // how many times ShellPoint called
00816   mySumDist += sqrt( sqDistance );
00817   MESSAGE ( " ------ SOLUTION: ( "<<solution.X()<<" "<<solution.Y()<<" "<<solution.Z()<<" )"<< std::endl
00818          << " ------ DIST : " << sqrt( sqDistance ) << "\t Tol=" << myTolerance << "\t Nb LOOPS=" << nbLoops << std::endl
00819          << " ------ NB IT: " << myNbIterations << ",  SUM DIST: " << mySumDist );
00820 #endif
00821 
00822   theParams = solution;
00823 
00824   if ( myFaceIndex > 0 )
00825     theParams.SetCoord( myFaceIndex, myFaceParam );
00826 
00827   return true;
00828 }
00829 
00830 //=======================================================================
00831 //function : VertexParameters
00832 //purpose  : return parameters of a vertex given by TShapeID
00833 //=======================================================================
00834 
00835 bool SMESH_Block::VertexParameters(const int theVertexID, gp_XYZ& theParams)
00836 {
00837   switch ( theVertexID ) {
00838   case ID_V000: theParams.SetCoord(0., 0., 0.); return true;
00839   case ID_V100: theParams.SetCoord(1., 0., 0.); return true;
00840   case ID_V110: theParams.SetCoord(1., 1., 0.); return true;
00841   case ID_V010: theParams.SetCoord(0., 1., 0.); return true;
00842   default:;
00843   }
00844   return false;
00845 }
00846 
00847 //=======================================================================
00848 //function : EdgeParameters
00849 //purpose  : return parameters of a point given by theU on edge
00850 //=======================================================================
00851 
00852 bool SMESH_Block::EdgeParameters(const int theEdgeID, const double theU, gp_XYZ& theParams)
00853 {
00854   if ( IsEdgeID( theEdgeID )) {
00855     vector< int > vertexVec;
00856     GetEdgeVertexIDs( theEdgeID, vertexVec );
00857     VertexParameters( vertexVec[0], theParams );
00858     TEdge& e = myEdge[ theEdgeID - ID_Ex00 ];
00859     double param = ( theU - e.EndParam(0) ) / ( e.EndParam(1) - e.EndParam(0) );
00860     theParams.SetCoord( e.CoordInd(), param );
00861     return true;
00862   }
00863   return false;
00864 }
00865 
00866 //=======================================================================
00867 //function : DumpShapeID
00868 //purpose  : debug an id of a block sub-shape
00869 //=======================================================================
00870 
00871 #define CASEDUMP(id,strm) case id: strm << #id; break;
00872 
00873 ostream& SMESH_Block::DumpShapeID (const int id, ostream& stream)
00874 {
00875   switch ( id ) {
00876   CASEDUMP( ID_V000, stream );
00877   CASEDUMP( ID_V100, stream );
00878   CASEDUMP( ID_V010, stream );
00879   CASEDUMP( ID_V110, stream );
00880   CASEDUMP( ID_V001, stream );
00881   CASEDUMP( ID_V101, stream );
00882   CASEDUMP( ID_V011, stream );
00883   CASEDUMP( ID_V111, stream );
00884   CASEDUMP( ID_Ex00, stream );
00885   CASEDUMP( ID_Ex10, stream );
00886   CASEDUMP( ID_Ex01, stream );
00887   CASEDUMP( ID_Ex11, stream );
00888   CASEDUMP( ID_E0y0, stream );
00889   CASEDUMP( ID_E1y0, stream );
00890   CASEDUMP( ID_E0y1, stream );
00891   CASEDUMP( ID_E1y1, stream );
00892   CASEDUMP( ID_E00z, stream );
00893   CASEDUMP( ID_E10z, stream );
00894   CASEDUMP( ID_E01z, stream );
00895   CASEDUMP( ID_E11z, stream );
00896   CASEDUMP( ID_Fxy0, stream );
00897   CASEDUMP( ID_Fxy1, stream );
00898   CASEDUMP( ID_Fx0z, stream );
00899   CASEDUMP( ID_Fx1z, stream );
00900   CASEDUMP( ID_F0yz, stream );
00901   CASEDUMP( ID_F1yz, stream );
00902   CASEDUMP( ID_Shell, stream );
00903   default: stream << "ID_INVALID";
00904   }
00905   return stream;
00906 }
00907 
00908 //=======================================================================
00909 //function : GetShapeIDByParams
00910 //purpose  : define an id of the block sub-shape by normlized point coord
00911 //=======================================================================
00912 
00913 int SMESH_Block::GetShapeIDByParams ( const gp_XYZ& theCoord )
00914 {
00915   //   id ( 0 - 26 ) computation:
00916 
00917   //   vertex     ( 0 - 7 )  : id = 1*x + 2*y + 4*z
00918 
00919   //   edge || X  ( 8 - 11 ) : id = 8   + 1*y + 2*z
00920   //   edge || Y  ( 12 - 15 ): id = 1*x + 12  + 2*z
00921   //   edge || Z  ( 16 - 19 ): id = 1*x + 2*y + 16 
00922 
00923   //   face || XY ( 20 - 21 ): id = 8   + 12  + 1*z - 0
00924   //   face || XZ ( 22 - 23 ): id = 8   + 1*y + 16  - 2
00925   //   face || YZ ( 24 - 25 ): id = 1*x + 12  + 16  - 4
00926 
00927   static int iAddBnd[]    = { 1, 2, 4 };
00928   static int iAddNotBnd[] = { 8, 12, 16 };
00929   static int iFaceSubst[] = { 0, 2, 4 };
00930 
00931   int id = 0;
00932   int iOnBoundary = 0;
00933   for ( int iCoord = 0; iCoord < 3; iCoord++ )
00934   {
00935     double val = theCoord.Coord( iCoord + 1 );
00936     if ( val == 0.0 )
00937       iOnBoundary++;
00938     else if ( val == 1.0 )
00939       id += iAddBnd[ iOnBoundary++ ];
00940     else
00941       id += iAddNotBnd[ iCoord ];
00942   }
00943   if ( iOnBoundary == 1 ) // face
00944     id -= iFaceSubst[ (id - 20) / 4 ];
00945   else if ( iOnBoundary == 0 ) // shell
00946     id = 26;
00947 
00948   if ( id > 26 || id < 0 ) {
00949     MESSAGE( "GetShapeIDByParams() = " << id
00950             <<" "<< theCoord.X() <<" "<< theCoord.Y() <<" "<< theCoord.Z() );
00951   }
00952 
00953   return id + 1; // shape ids start at 1
00954 }
00955 
00956 //================================================================================
00972 //================================================================================
00973 
00974 int SMESH_Block::GetOrderedEdges (const TopoDS_Face&   theFace,
00975                                   TopoDS_Vertex        theFirstVertex,
00976                                   list< TopoDS_Edge >& theEdges,
00977                                   list< int >  &       theNbEdgesInWires,
00978                                   const bool           theShapeAnalysisAlgo)
00979 {
00980   // put wires in a list, so that an outer wire comes first
00981   list<TopoDS_Wire> aWireList;
00982   TopoDS_Wire anOuterWire =
00983     theShapeAnalysisAlgo ? ShapeAnalysis::OuterWire( theFace ) : BRepTools::OuterWire( theFace );
00984   for ( TopoDS_Iterator wIt (theFace); wIt.More(); wIt.Next() )
00985     if ( wIt.Value().ShapeType() == TopAbs_WIRE ) // it can be internal vertex!
00986     {
00987       if ( !anOuterWire.IsSame( wIt.Value() ))
00988         aWireList.push_back( TopoDS::Wire( wIt.Value() ));
00989       else
00990         aWireList.push_front( TopoDS::Wire( wIt.Value() ));
00991     }
00992 
00993   // loop on edges of wires
00994   theNbEdgesInWires.clear();
00995   list<TopoDS_Wire>::iterator wlIt = aWireList.begin();
00996   for ( ; wlIt != aWireList.end(); wlIt++ )
00997   {
00998     int iE;
00999     BRepTools_WireExplorer wExp( *wlIt, theFace );
01000     for ( iE = 0; wExp.More(); wExp.Next(), iE++ )
01001     {
01002       TopoDS_Edge edge = wExp.Current();
01003       // commented for issue 0020557, other related ones: 0020526, PAL19080
01004       // edge = TopoDS::Edge( edge.Oriented( wExp.Orientation() ));
01005       theEdges.push_back( edge );
01006     }
01007     if ( iE == 0 ) // wExp returns nothing if e.g. the wire contains one internal edge
01008     { // Issue 0020676
01009       for ( TopoDS_Iterator e( *wlIt ); e.More(); e.Next(), ++iE )
01010         theEdges.push_back( TopoDS::Edge( e.Value() ));
01011     }
01012     theNbEdgesInWires.push_back( iE );
01013     iE = 0;
01014     if ( wlIt == aWireList.begin() && theEdges.size() > 1 ) { // the outer wire
01015       // orient closed edges
01016       list< TopoDS_Edge >::iterator eIt, eIt2;
01017       for ( eIt = theEdges.begin(); eIt != theEdges.end(); eIt++ )
01018       {
01019         TopoDS_Edge& edge = *eIt;
01020         if ( TopExp::FirstVertex( edge ).IsSame( TopExp::LastVertex( edge ) ))
01021         {
01022           eIt2 = eIt;
01023           bool isNext = ( eIt2 == theEdges.begin() );
01024           TopoDS_Edge edge2 = isNext ? *(++eIt2) : *(--eIt2);
01025           double f1,l1,f2,l2;
01026           Handle(Geom2d_Curve) c1 = BRep_Tool::CurveOnSurface( edge, theFace, f1,l1 );
01027           Handle(Geom2d_Curve) c2 = BRep_Tool::CurveOnSurface( edge2, theFace, f2,l2 );
01028           gp_Pnt2d pf = c1->Value( edge.Orientation() == TopAbs_FORWARD ? f1 : l1 );
01029           gp_Pnt2d pl = c1->Value( edge.Orientation() == TopAbs_FORWARD ? l1 : f1 );
01030           bool isFirst = ( edge2.Orientation() == TopAbs_FORWARD ? isNext : !isNext );
01031           gp_Pnt2d p2 = c2->Value( isFirst ? f2 : l2 );
01032           isFirst = ( p2.SquareDistance( pf ) < p2.SquareDistance( pl ));
01033           if ( isNext ? isFirst : !isFirst )
01034             edge.Reverse();
01035           // to make a seam go first
01036           if ( theFirstVertex.IsNull() )
01037             theFirstVertex = TopExp::FirstVertex( edge, true );
01038         }
01039       }
01040       // rotate theEdges until it begins from theFirstVertex
01041       if ( ! theFirstVertex.IsNull() ) {
01042         TopoDS_Vertex vv[2];
01043         TopExp::Vertices( theEdges.front(), vv[0], vv[1], true );
01044         // on closed face, make seam edge the first in the list
01045         while ( !vv[0].IsSame( theFirstVertex ) || vv[0].IsSame( vv[1] ))
01046         {
01047           theEdges.splice(theEdges.end(), theEdges,
01048                           theEdges.begin(), ++theEdges.begin());
01049           TopExp::Vertices( theEdges.front(), vv[0], vv[1], true );
01050           if ( iE++ > theNbEdgesInWires.back() ) {
01051 #ifdef _DEBUG_
01052             gp_Pnt p = BRep_Tool::Pnt( theFirstVertex );
01053             MESSAGE ( " : Warning : vertex "<< theFirstVertex.TShape().operator->()
01054                    << " ( " << p.X() << " " << p.Y() << " " << p.Z() << " )" 
01055                    << " not found in outer wire of face "<< theFace.TShape().operator->()
01056                    << " with vertices: " );
01057             wExp.Init( *wlIt, theFace );
01058             for ( int i = 0; wExp.More(); wExp.Next(), i++ )
01059             {
01060               TopoDS_Edge edge = wExp.Current();
01061               edge = TopoDS::Edge( edge.Oriented( wExp.Orientation() ));
01062               TopoDS_Vertex v = TopExp::FirstVertex( edge, true );
01063               gp_Pnt p = BRep_Tool::Pnt( v );
01064               MESSAGE_ADD ( i << " " << v.TShape().operator->() << " "
01065                             << p.X() << " " << p.Y() << " " << p.Z() << " " << std::endl );
01066             }
01067 #endif
01068             break; // break infinite loop
01069           }
01070         }
01071       }
01072     } // end outer wire
01073   }
01074 
01075   return aWireList.size();
01076 }
01077 //================================================================================
01081 //================================================================================
01082 
01083 void SMESH_Block::init()
01084 {
01085   myNbIterations = 0;
01086   mySumDist = 0;
01087   myGridComputed = false;
01088 }
01089 
01090 //=======================================================================
01091 //function : LoadMeshBlock
01092 //purpose  : prepare to work with theVolume
01093 //=======================================================================
01094 
01095 #define gpXYZ(n) gp_XYZ(n->X(),n->Y(),n->Z())
01096 
01097 bool SMESH_Block::LoadMeshBlock(const SMDS_MeshVolume*        theVolume,
01098                                 const int                     theNode000Index,
01099                                 const int                     theNode001Index,
01100                                 vector<const SMDS_MeshNode*>& theOrderedNodes)
01101 {
01102   MESSAGE(" ::LoadMeshBlock()");
01103   init();
01104 
01105   SMDS_VolumeTool vTool;
01106   if (!vTool.Set( theVolume ) || vTool.NbNodes() != 8 ||
01107       !vTool.IsLinked( theNode000Index, theNode001Index )) {
01108     MESSAGE(" Bad arguments ");
01109     return false;
01110   }
01111   vTool.SetExternalNormal();
01112   // In terms of indices used for access to nodes and faces in SMDS_VolumeTool:
01113   int V000, V100, V010, V110, V001, V101, V011, V111; // 8 vertices
01114   int Fxy0, Fxy1; // bottom and top faces
01115   // vertices of faces
01116   vector<int> vFxy0, vFxy1;
01117 
01118   V000 = theNode000Index;
01119   V001 = theNode001Index;
01120 
01121   // get faces sharing V000 and V001
01122   list<int> fV000, fV001;
01123   int i, iF, iE, iN;
01124   for ( iF = 0; iF < vTool.NbFaces(); ++iF ) {
01125     const int* nid = vTool.GetFaceNodesIndices( iF );
01126     for ( iN = 0; iN < 4; ++iN )
01127       if ( nid[ iN ] == V000 ) {
01128         fV000.push_back( iF );
01129       } else if ( nid[ iN ] == V001 ) {
01130         fV001.push_back( iF );
01131       }
01132   }
01133 
01134   // find the bottom (Fxy0), the top (Fxy1) faces
01135   list<int>::iterator fIt1, fIt2, Fxy0Pos;
01136   for ( fIt1 = fV000.begin(); fIt1 != fV000.end(); fIt1++) {
01137     fIt2 = std::find( fV001.begin(), fV001.end(), *fIt1 );
01138     if ( fIt2 != fV001.end() ) { // *fIt1 is in the both lists
01139       fV001.erase( fIt2 ); // erase Fx0z or F0yz from fV001
01140     } else { // *fIt1 is in fV000 only
01141       Fxy0Pos = fIt1; // points to Fxy0
01142     }
01143   }
01144   Fxy0 = *Fxy0Pos;
01145   Fxy1 = fV001.front();
01146   const SMDS_MeshNode** nn = vTool.GetNodes();
01147 
01148   // find bottom veritices, their order is that a face normal is external
01149   vFxy0.resize(4);
01150   const int* nid = vTool.GetFaceNodesIndices( Fxy0 );
01151   for ( i = 0; i < 4; ++i )
01152     if ( nid[ i ] == V000 )
01153       break;
01154   for ( iN = 0; iN < 4; ++iN, ++i ) {
01155     if ( i == 4 ) i = 0;
01156     vFxy0[ iN ] = nid[ i ];
01157   }
01158   // find top veritices, their order is that a face normal is external
01159   vFxy1.resize(4);
01160   nid = vTool.GetFaceNodesIndices( Fxy1 );
01161   for ( i = 0; i < 4; ++i )
01162     if ( nid[ i ] == V001 )
01163       break;
01164   for ( iN = 0; iN < 4; ++iN, ++i ) {
01165     if ( i == 4 ) i = 0;
01166     vFxy1[ iN ] = nid[ i ];
01167   }
01168   // find indices of the rest veritices 
01169   V100 = vFxy0[3];
01170   V010 = vFxy0[1];
01171   V110 = vFxy0[2];
01172   V101 = vFxy1[1];
01173   V011 = vFxy1[3];
01174   V111 = vFxy1[2];
01175 
01176   // set points coordinates
01177   myPnt[ ID_V000 - 1 ] = gpXYZ( nn[ V000 ] );
01178   myPnt[ ID_V100 - 1 ] = gpXYZ( nn[ V100 ] );
01179   myPnt[ ID_V010 - 1 ] = gpXYZ( nn[ V010 ] );
01180   myPnt[ ID_V110 - 1 ] = gpXYZ( nn[ V110 ] );
01181   myPnt[ ID_V001 - 1 ] = gpXYZ( nn[ V001 ] );
01182   myPnt[ ID_V101 - 1 ] = gpXYZ( nn[ V101 ] );
01183   myPnt[ ID_V011 - 1 ] = gpXYZ( nn[ V011 ] );
01184   myPnt[ ID_V111 - 1 ] = gpXYZ( nn[ V111 ] );
01185 
01186   // fill theOrderedNodes
01187   theOrderedNodes.resize( 8 );
01188   theOrderedNodes[ 0 ] = nn[ V000 ];
01189   theOrderedNodes[ 1 ] = nn[ V100 ];
01190   theOrderedNodes[ 2 ] = nn[ V010 ];
01191   theOrderedNodes[ 3 ] = nn[ V110 ];
01192   theOrderedNodes[ 4 ] = nn[ V001 ];
01193   theOrderedNodes[ 5 ] = nn[ V101 ];
01194   theOrderedNodes[ 6 ] = nn[ V011 ];
01195   theOrderedNodes[ 7 ] = nn[ V111 ];
01196   
01197   // fill edges
01198   vector< int > vertexVec;
01199   for ( iE = 0; iE < NbEdges(); ++iE ) {
01200     GetEdgeVertexIDs(( iE + ID_FirstE ), vertexVec );
01201     myEdge[ iE ].Set(( iE + ID_FirstE ),
01202                      myPnt[ vertexVec[0] - 1 ],
01203                      myPnt[ vertexVec[1] - 1 ]);
01204   }
01205 
01206   // fill faces' corners
01207   for ( iF = ID_Fxy0; iF < ID_Shell; ++iF )
01208   {
01209     TFace& tFace = myFace[ iF - ID_FirstF ];
01210     vector< int > edgeIdVec(4, -1);
01211     GetFaceEdgesIDs( iF, edgeIdVec );
01212     tFace.Set( iF, myEdge[ edgeIdVec [ 0 ] - ID_Ex00], myEdge[ edgeIdVec [ 1 ] - ID_Ex00]);
01213   }
01214 
01215   return true;
01216 }
01217 
01218 //=======================================================================
01219 //function : LoadBlockShapes
01220 //purpose  : Initialize block geometry with theShell,
01221 //           add sub-shapes of theBlock to theShapeIDMap so that they get
01222 //           IDs acoording to enum TShapeID
01223 //=======================================================================
01224 
01225 bool SMESH_Block::LoadBlockShapes(const TopoDS_Shell&         theShell,
01226                                   const TopoDS_Vertex&        theVertex000,
01227                                   const TopoDS_Vertex&        theVertex001,
01228                                   TopTools_IndexedMapOfOrientedShape& theShapeIDMap )
01229 {
01230   MESSAGE(" ::LoadBlockShapes()");
01231   return ( FindBlockShapes( theShell, theVertex000, theVertex001, theShapeIDMap ) &&
01232            LoadBlockShapes( theShapeIDMap ));
01233 }
01234 
01235 //=======================================================================
01236 //function : LoadBlockShapes
01237 //purpose  : add sub-shapes of theBlock to theShapeIDMap so that they get
01238 //           IDs acoording to enum TShapeID
01239 //=======================================================================
01240 
01241 bool SMESH_Block::FindBlockShapes(const TopoDS_Shell&         theShell,
01242                                   const TopoDS_Vertex&        theVertex000,
01243                                   const TopoDS_Vertex&        theVertex001,
01244                                   TopTools_IndexedMapOfOrientedShape& theShapeIDMap )
01245 {
01246   MESSAGE(" ::FindBlockShapes()");
01247 
01248   // 8 vertices
01249   TopoDS_Shape V000, V100, V010, V110, V001, V101, V011, V111;
01250   // 12 edges
01251   TopoDS_Shape Ex00, Ex10, Ex01, Ex11;
01252   TopoDS_Shape E0y0, E1y0, E0y1, E1y1;
01253   TopoDS_Shape E00z, E10z, E01z, E11z;
01254   // 6 faces
01255   TopoDS_Shape Fxy0, Fx0z, F0yz, Fxy1, Fx1z, F1yz;
01256 
01257   // nb of faces bound to a vertex in TopTools_IndexedDataMapOfShapeListOfShape
01258   // filled by TopExp::MapShapesAndAncestors()
01259   const int NB_FACES_BY_VERTEX = 6;
01260 
01261   TopTools_IndexedDataMapOfShapeListOfShape vfMap;
01262   TopExp::MapShapesAndAncestors( theShell, TopAbs_VERTEX, TopAbs_FACE, vfMap );
01263   if ( vfMap.Extent() != 8 ) {
01264     MESSAGE(" Wrong nb of vertices in the block: " << vfMap.Extent() );
01265     return false;
01266   }
01267 
01268   V000 = theVertex000;
01269   V001 = theVertex001;
01270 
01271   if ( V000.IsNull() ) {
01272     // find vertex 000 - the one with smallest coordinates
01273     double minVal = DBL_MAX, minX, val;
01274     for ( int i = 1; i <= 8; i++ ) {
01275       const TopoDS_Vertex& v = TopoDS::Vertex( vfMap.FindKey( i ));
01276       gp_Pnt P = BRep_Tool::Pnt( v );
01277       val = P.X() + P.Y() + P.Z();
01278       if ( val < minVal || ( val == minVal && P.X() < minX )) {
01279         V000 = v;
01280         minVal = val;
01281         minX = P.X();
01282       }
01283     }
01284     // find vertex 001 - the one on the most vertical edge passing through V000
01285     TopTools_IndexedDataMapOfShapeListOfShape veMap;
01286     TopExp::MapShapesAndAncestors( theShell, TopAbs_VERTEX, TopAbs_EDGE, veMap );
01287     gp_Vec dir001 = gp::DZ();
01288     gp_Pnt p000 = BRep_Tool::Pnt( TopoDS::Vertex( V000 ));
01289     double maxVal = -DBL_MAX;
01290     TopTools_ListIteratorOfListOfShape eIt ( veMap.FindFromKey( V000 ));
01291     for (  ; eIt.More(); eIt.Next() ) {
01292       const TopoDS_Edge& e = TopoDS::Edge( eIt.Value() );
01293       TopoDS_Vertex v = TopExp::FirstVertex( e );
01294       if ( v.IsSame( V000 ))
01295         v = TopExp::LastVertex( e );
01296       val = dir001 * gp_Vec( p000, BRep_Tool::Pnt( v )).Normalized();
01297       if ( val > maxVal ) {
01298         V001 = v;
01299         maxVal = val;
01300       }
01301     }
01302   }
01303 
01304   // find the bottom (Fxy0), Fx0z and F0yz faces
01305 
01306   const TopTools_ListOfShape& f000List = vfMap.FindFromKey( V000 );
01307   const TopTools_ListOfShape& f001List = vfMap.FindFromKey( V001 );
01308   if (f000List.Extent() != NB_FACES_BY_VERTEX ||
01309       f001List.Extent() != NB_FACES_BY_VERTEX ) {
01310     MESSAGE(" LoadBlockShapes() " << f000List.Extent() << " " << f001List.Extent());
01311     return false;
01312   }
01313   TopTools_ListIteratorOfListOfShape f001It, f000It ( f000List );
01314   int i, j, iFound1, iFound2;
01315   for ( j = 0; f000It.More(); f000It.Next(), j++ )
01316   {
01317     if ( NB_FACES_BY_VERTEX == 6 && j % 2 ) continue; // each face encounters twice
01318     const TopoDS_Shape& F = f000It.Value();
01319     for ( i = 0, f001It.Initialize( f001List ); f001It.More(); f001It.Next(), i++ ) {
01320       if ( NB_FACES_BY_VERTEX == 6 && i % 2 ) continue; // each face encounters twice
01321       if ( F.IsSame( f001It.Value() ))
01322         break;
01323     }
01324     if ( f001It.More() ) // Fx0z or F0yz found
01325       if ( Fx0z.IsNull() ) {
01326         Fx0z = F;
01327         iFound1 = i;
01328       } else {
01329         F0yz = F;
01330         iFound2 = i;
01331       }
01332     else // F is the bottom face
01333       Fxy0 = F;
01334   }
01335   if ( Fxy0.IsNull() || Fx0z.IsNull() || F0yz.IsNull() ) {
01336     MESSAGE( Fxy0.IsNull() <<" "<< Fx0z.IsNull() <<" "<< F0yz.IsNull() );
01337     return false;
01338   }
01339 
01340   // choose the top face (Fxy1)
01341   for ( i = 0, f001It.Initialize( f001List ); f001It.More(); f001It.Next(), i++ ) {
01342     if ( NB_FACES_BY_VERTEX == 6 && i % 2 ) continue; // each face encounters twice
01343     if ( i != iFound1 && i != iFound2 )
01344       break;
01345   }
01346   Fxy1 = f001It.Value();
01347   if ( Fxy1.IsNull() ) {
01348     MESSAGE(" LoadBlockShapes() error ");
01349     return false;
01350   }
01351 
01352   // find bottom edges and veritices
01353   list< TopoDS_Edge > eList;
01354   list< int >         nbVertexInWires;
01355   GetOrderedEdges( TopoDS::Face( Fxy0 ), TopoDS::Vertex( V000 ), eList, nbVertexInWires );
01356   if ( nbVertexInWires.size() != 1 || nbVertexInWires.front() != 4 ) {
01357     MESSAGE(" LoadBlockShapes() error ");
01358     return false;
01359   }
01360   list< TopoDS_Edge >::iterator elIt = eList.begin();
01361   for ( i = 0; elIt != eList.end(); elIt++, i++ )
01362     switch ( i ) {
01363     case 0: E0y0 = *elIt; V010 = TopExp::LastVertex( *elIt, true ); break;
01364     case 1: Ex10 = *elIt; V110 = TopExp::LastVertex( *elIt, true ); break;
01365     case 2: E1y0 = *elIt; V100 = TopExp::LastVertex( *elIt, true ); break;
01366     case 3: Ex00 = *elIt; break;
01367     default:;
01368     }
01369   if ( i != 4 || E0y0.IsNull() || Ex10.IsNull() || E1y0.IsNull() || Ex00.IsNull() ) {
01370     MESSAGE(" LoadBlockShapes() error, eList.size()=" << eList.size());
01371     return false;
01372   }
01373 
01374 
01375   // find top edges and veritices
01376   eList.clear();
01377   GetOrderedEdges( TopoDS::Face( Fxy1 ), TopoDS::Vertex( V001 ), eList, nbVertexInWires );
01378   if ( nbVertexInWires.size() != 1 || nbVertexInWires.front() != 4 ) {
01379     MESSAGE(" LoadBlockShapes() error ");
01380     return false;
01381   }
01382   for ( i = 0, elIt = eList.begin(); elIt != eList.end(); elIt++, i++ )
01383     switch ( i ) {
01384     case 0: Ex01 = *elIt; V101 = TopExp::LastVertex( *elIt, true ); break;
01385     case 1: E1y1 = *elIt; V111 = TopExp::LastVertex( *elIt, true ); break;
01386     case 2: Ex11 = *elIt; V011 = TopExp::LastVertex( *elIt, true ); break;
01387     case 3: E0y1 = *elIt; break;
01388     default:;
01389     }
01390   if ( i != 4 || Ex01.IsNull() || E1y1.IsNull() || Ex11.IsNull() || E0y1.IsNull() ) {
01391     MESSAGE(" LoadBlockShapes() error, eList.size()=" << eList.size());
01392     return false;
01393   }
01394 
01395   // swap Fx0z and F0yz if necessary
01396   TopExp_Explorer exp( Fx0z, TopAbs_VERTEX );
01397   for ( ; exp.More(); exp.Next() ) // Fx0z shares V101 and V100
01398     if ( V101.IsSame( exp.Current() ) || V100.IsSame( exp.Current() ))
01399       break; // V101 or V100 found
01400   if ( !exp.More() ) { // not found
01401     std::swap( Fx0z, F0yz);
01402   }
01403 
01404   // find Fx1z and F1yz faces
01405   const TopTools_ListOfShape& f111List = vfMap.FindFromKey( V111 );
01406   const TopTools_ListOfShape& f110List = vfMap.FindFromKey( V110 );
01407   if (f111List.Extent() != NB_FACES_BY_VERTEX ||
01408       f110List.Extent() != NB_FACES_BY_VERTEX ) {
01409     MESSAGE(" LoadBlockShapes() " << f111List.Extent() << " " << f110List.Extent());
01410     return false;
01411   }
01412   TopTools_ListIteratorOfListOfShape f111It, f110It ( f110List);
01413   for ( j = 0 ; f110It.More(); f110It.Next(), j++ ) {
01414     if ( NB_FACES_BY_VERTEX == 6 && j % 2 ) continue; // each face encounters twice
01415     const TopoDS_Shape& F = f110It.Value();
01416     for ( i = 0, f111It.Initialize( f111List ); f111It.More(); f111It.Next(), i++ ) {
01417       if ( NB_FACES_BY_VERTEX == 6 && i % 2 ) continue; // each face encounters twice
01418       if ( F.IsSame( f111It.Value() )) { // Fx1z or F1yz found
01419         if ( Fx1z.IsNull() )
01420           Fx1z = F;
01421         else
01422           F1yz = F;
01423       }
01424     }
01425   }
01426   if ( Fx1z.IsNull() || F1yz.IsNull() ) {
01427     MESSAGE(" LoadBlockShapes() error ");
01428     return false;
01429   }
01430 
01431   // swap Fx1z and F1yz if necessary
01432   for ( exp.Init( Fx1z, TopAbs_VERTEX ); exp.More(); exp.Next() )
01433     if ( V010.IsSame( exp.Current() ) || V011.IsSame( exp.Current() ))
01434       break;
01435   if ( !exp.More() ) {
01436     std::swap( Fx1z, F1yz);
01437   }
01438 
01439   // find vertical edges
01440   for ( exp.Init( Fx0z, TopAbs_EDGE ); exp.More(); exp.Next() ) {
01441     const TopoDS_Edge& edge = TopoDS::Edge( exp.Current() );
01442     const TopoDS_Shape& vFirst = TopExp::FirstVertex( edge, true );
01443     if ( vFirst.IsSame( V001 ))
01444       E00z = edge;
01445     else if ( vFirst.IsSame( V100 ))
01446       E10z = edge;
01447   }
01448   if ( E00z.IsNull() || E10z.IsNull() ) {
01449     MESSAGE(" LoadBlockShapes() error ");
01450     return false;
01451   }
01452   for ( exp.Init( Fx1z, TopAbs_EDGE ); exp.More(); exp.Next() ) {
01453     const TopoDS_Edge& edge = TopoDS::Edge( exp.Current() );
01454     const TopoDS_Shape& vFirst = TopExp::FirstVertex( edge, true );
01455     if ( vFirst.IsSame( V111 ))
01456       E11z = edge;
01457     else if ( vFirst.IsSame( V010 ))
01458       E01z = edge;
01459   }
01460   if ( E01z.IsNull() || E11z.IsNull() ) {
01461     MESSAGE(" LoadBlockShapes() error ");
01462     return false;
01463   }
01464 
01465   // load shapes in theShapeIDMap
01466 
01467   theShapeIDMap.Clear();
01468   
01469   theShapeIDMap.Add(V000.Oriented( TopAbs_FORWARD ));
01470   theShapeIDMap.Add(V100.Oriented( TopAbs_FORWARD ));
01471   theShapeIDMap.Add(V010.Oriented( TopAbs_FORWARD ));
01472   theShapeIDMap.Add(V110.Oriented( TopAbs_FORWARD ));
01473   theShapeIDMap.Add(V001.Oriented( TopAbs_FORWARD ));
01474   theShapeIDMap.Add(V101.Oriented( TopAbs_FORWARD ));
01475   theShapeIDMap.Add(V011.Oriented( TopAbs_FORWARD ));
01476   theShapeIDMap.Add(V111.Oriented( TopAbs_FORWARD ));
01477 
01478   theShapeIDMap.Add(Ex00);
01479   theShapeIDMap.Add(Ex10);
01480   theShapeIDMap.Add(Ex01);
01481   theShapeIDMap.Add(Ex11);
01482 
01483   theShapeIDMap.Add(E0y0);
01484   theShapeIDMap.Add(E1y0);
01485   theShapeIDMap.Add(E0y1);
01486   theShapeIDMap.Add(E1y1);
01487 
01488   theShapeIDMap.Add(E00z);
01489   theShapeIDMap.Add(E10z);
01490   theShapeIDMap.Add(E01z);
01491   theShapeIDMap.Add(E11z);
01492 
01493   theShapeIDMap.Add(Fxy0);
01494   theShapeIDMap.Add(Fxy1);
01495   theShapeIDMap.Add(Fx0z);
01496   theShapeIDMap.Add(Fx1z);
01497   theShapeIDMap.Add(F0yz);
01498   theShapeIDMap.Add(F1yz);
01499   
01500   theShapeIDMap.Add(theShell);
01501 
01502   return true;
01503 }
01504 
01505 //================================================================================
01511 //================================================================================
01512 
01513 bool SMESH_Block::LoadBlockShapes(const TopTools_IndexedMapOfOrientedShape& theShapeIDMap)
01514 {
01515   init();
01516 
01517   // store shapes geometry
01518   for ( int shapeID = 1; shapeID < theShapeIDMap.Extent(); shapeID++ )
01519   {
01520     const TopoDS_Shape& S = theShapeIDMap( shapeID );
01521     switch ( S.ShapeType() )
01522     {
01523     case TopAbs_VERTEX: {
01524 
01525       if ( !IsVertexID( ID_V111 )) return false;
01526       myPnt[ shapeID - ID_V000 ] = BRep_Tool::Pnt( TopoDS::Vertex( S )).XYZ();
01527       break;
01528     }
01529     case TopAbs_EDGE: {
01530 
01531       if ( !IsEdgeID( shapeID )) return false;
01532       const TopoDS_Edge& edge = TopoDS::Edge( S );
01533       TEdge& tEdge = myEdge[ shapeID - ID_FirstE ];
01534       tEdge.Set( shapeID,
01535                  new BRepAdaptor_Curve( edge ),
01536                  IsForwardEdge( edge, theShapeIDMap ));
01537       break;
01538     }
01539     case TopAbs_FACE: {
01540 
01541       if ( !LoadFace( TopoDS::Face( S ), shapeID, theShapeIDMap ))
01542         return false;
01543       break;
01544     }
01545     default: break;
01546     }
01547   } // loop on shapes in theShapeIDMap
01548 
01549   return true;
01550 }
01551 
01552 //================================================================================
01563 //================================================================================
01564 
01565 bool SMESH_Block::LoadFace(const TopoDS_Face& theFace,
01566                            const int          theFaceID,
01567                            const TopTools_IndexedMapOfOrientedShape& theShapeIDMap)
01568 {
01569   if ( !IsFaceID( theFaceID ) ) return false;
01570   // pcurves
01571   Adaptor2d_Curve2d* c2d[4];
01572   bool isForward[4];
01573   vector< int > edgeIdVec;
01574   GetFaceEdgesIDs( theFaceID, edgeIdVec );
01575   for ( int iE = 0; iE < edgeIdVec.size(); iE++ ) // loop on 4 edges
01576   {
01577     if ( edgeIdVec[ iE ] > theShapeIDMap.Extent() )
01578       return false;
01579     const TopoDS_Edge& edge = TopoDS::Edge( theShapeIDMap( edgeIdVec[ iE ]));
01580     c2d[ iE ] = new BRepAdaptor_Curve2d( edge, theFace );
01581     isForward[ iE ] = IsForwardEdge( edge, theShapeIDMap );
01582   }
01583   TFace& tFace = myFace[ theFaceID - ID_FirstF ];
01584   tFace.Set( theFaceID, new BRepAdaptor_Surface( theFace ), c2d, isForward );
01585   return true;
01586 }
01587 
01588 //================================================================================
01595 //================================================================================
01596 
01597 bool SMESH_Block::Insert(const TopoDS_Shape& theShape,
01598                          const int           theShapeID,
01599                          TopTools_IndexedMapOfOrientedShape& theShapeIDMap)
01600 {
01601   if ( !theShape.IsNull() && theShapeID > 0 )
01602   {
01603     if ( theShapeIDMap.Contains( theShape ))
01604       return ( theShapeIDMap.FindIndex( theShape ) == theShapeID );
01605 
01606     if ( theShapeID <= theShapeIDMap.Extent() ) {
01607         theShapeIDMap.Substitute( theShapeID, theShape );
01608     }
01609     else {
01610       while ( theShapeIDMap.Extent() < theShapeID - 1 ) {
01611         TopoDS_Compound comp;
01612         BRep_Builder().MakeCompound( comp );
01613         theShapeIDMap.Add( comp );
01614       }
01615       theShapeIDMap.Add( theShape );
01616     }
01617     return true;
01618   }
01619   return false;
01620 }
01621 
01622 //=======================================================================
01623 //function : GetFaceEdgesIDs
01624 //purpose  : return edges IDs in the order u0, u1, 0v, 1v
01625 //           u0 means "|| u, v == 0"
01626 //=======================================================================
01627 
01628 void SMESH_Block::GetFaceEdgesIDs (const int faceID, vector< int >& edgeVec )
01629 {
01630   edgeVec.resize( 4 );
01631   switch ( faceID ) {
01632   case ID_Fxy0:
01633     edgeVec[ 0 ] = ID_Ex00;
01634     edgeVec[ 1 ] = ID_Ex10;
01635     edgeVec[ 2 ] = ID_E0y0;
01636     edgeVec[ 3 ] = ID_E1y0;
01637     break;
01638   case ID_Fxy1:
01639     edgeVec[ 0 ] = ID_Ex01;
01640     edgeVec[ 1 ] = ID_Ex11;
01641     edgeVec[ 2 ] = ID_E0y1;
01642     edgeVec[ 3 ] = ID_E1y1;
01643     break;
01644   case ID_Fx0z:
01645     edgeVec[ 0 ] = ID_Ex00;
01646     edgeVec[ 1 ] = ID_Ex01;
01647     edgeVec[ 2 ] = ID_E00z;
01648     edgeVec[ 3 ] = ID_E10z;
01649     break;
01650   case ID_Fx1z:
01651     edgeVec[ 0 ] = ID_Ex10;
01652     edgeVec[ 1 ] = ID_Ex11;
01653     edgeVec[ 2 ] = ID_E01z;
01654     edgeVec[ 3 ] = ID_E11z;
01655     break;
01656   case ID_F0yz:
01657     edgeVec[ 0 ] = ID_E0y0;
01658     edgeVec[ 1 ] = ID_E0y1;
01659     edgeVec[ 2 ] = ID_E00z;
01660     edgeVec[ 3 ] = ID_E01z;
01661     break;
01662   case ID_F1yz:
01663     edgeVec[ 0 ] = ID_E1y0;
01664     edgeVec[ 1 ] = ID_E1y1;
01665     edgeVec[ 2 ] = ID_E10z;
01666     edgeVec[ 3 ] = ID_E11z;
01667     break;
01668   default:
01669     MESSAGE(" GetFaceEdgesIDs(), wrong face ID: " << faceID );
01670   }
01671 }
01672 
01673 //=======================================================================
01674 //function : GetEdgeVertexIDs
01675 //purpose  : return vertex IDs of an edge
01676 //=======================================================================
01677 
01678 void SMESH_Block::GetEdgeVertexIDs (const int edgeID, vector< int >& vertexVec )
01679 {
01680   vertexVec.resize( 2 );
01681   switch ( edgeID ) {
01682 
01683   case ID_Ex00:
01684     vertexVec[ 0 ] = ID_V000;
01685     vertexVec[ 1 ] = ID_V100;
01686     break;
01687   case ID_Ex10:
01688     vertexVec[ 0 ] = ID_V010;
01689     vertexVec[ 1 ] = ID_V110;
01690     break;
01691   case ID_Ex01:
01692     vertexVec[ 0 ] = ID_V001;
01693     vertexVec[ 1 ] = ID_V101;
01694     break;
01695   case ID_Ex11:
01696     vertexVec[ 0 ] = ID_V011;
01697     vertexVec[ 1 ] = ID_V111;
01698     break;
01699 
01700   case ID_E0y0:
01701     vertexVec[ 0 ] = ID_V000;
01702     vertexVec[ 1 ] = ID_V010;
01703     break;
01704   case ID_E1y0:
01705     vertexVec[ 0 ] = ID_V100;
01706     vertexVec[ 1 ] = ID_V110;
01707     break;
01708   case ID_E0y1:
01709     vertexVec[ 0 ] = ID_V001;
01710     vertexVec[ 1 ] = ID_V011;
01711     break;
01712   case ID_E1y1:
01713     vertexVec[ 0 ] = ID_V101;
01714     vertexVec[ 1 ] = ID_V111;
01715     break;
01716 
01717   case ID_E00z:
01718     vertexVec[ 0 ] = ID_V000;
01719     vertexVec[ 1 ] = ID_V001;
01720     break;
01721   case ID_E10z:
01722     vertexVec[ 0 ] = ID_V100;
01723     vertexVec[ 1 ] = ID_V101;
01724     break;
01725   case ID_E01z:
01726     vertexVec[ 0 ] = ID_V010;
01727     vertexVec[ 1 ] = ID_V011;
01728     break;
01729   case ID_E11z:
01730     vertexVec[ 0 ] = ID_V110;
01731     vertexVec[ 1 ] = ID_V111;
01732     break;
01733   default:
01734     vertexVec.resize(0);
01735     MESSAGE(" GetEdgeVertexIDs(), wrong edge ID: " << edgeID );
01736   }
01737 }
Copyright © 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE
Copyright © 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS