Back to index

salome-smesh  6.5.0
StdMeshers_SegmentLengthAroundVertex.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  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 //  SMESH SMESH : implementaion of SMESH idl descriptions
00024 //  File   : StdMeshers_SegmentLengthAroundVertex.cxx
00025 //  Module : SMESH
00026 //
00027 #include "StdMeshers_SegmentLengthAroundVertex.hxx"
00028 
00029 #include "SMESH_Mesh.hxx"
00030 #include "SMESH_Algo.hxx"
00031 #include "SMDS_MeshNode.hxx"
00032 #include "SMESHDS_Mesh.hxx"
00033 #include "SMESHDS_SubMesh.hxx"
00034 #include "SMESH_MesherHelper.hxx"
00035 
00036 #include "utilities.h"
00037 
00038 #include <BRepAdaptor_Curve.hxx>
00039 #include <GCPnts_AbscissaPoint.hxx>
00040 #include <TopTools_IndexedMapOfShape.hxx>
00041 #include <TopoDS.hxx>
00042 #include <TopoDS_Edge.hxx>
00043 
00044 using namespace std;
00045 
00046 //=============================================================================
00050 //=============================================================================
00051 
00052 StdMeshers_SegmentLengthAroundVertex::StdMeshers_SegmentLengthAroundVertex
00053                                        (int hypId, int studyId, SMESH_Gen * gen)
00054   :SMESH_Hypothesis(hypId, studyId, gen)
00055 {
00056   _length = 1.;
00057   _name = "SegmentLengthAroundVertex";
00058   _param_algo_dim = 0; // is used by StdMeshers_SegmentAroundVertex_0D
00059 }
00060 
00061 //=============================================================================
00065 //=============================================================================
00066 
00067 StdMeshers_SegmentLengthAroundVertex::~StdMeshers_SegmentLengthAroundVertex()
00068 {
00069 }
00070 
00071 //=============================================================================
00075 //=============================================================================
00076 
00077 void StdMeshers_SegmentLengthAroundVertex::SetLength(double length) throw(SALOME_Exception)
00078 {
00079   if (length <= 0)
00080     throw SALOME_Exception(LOCALIZED("length must be positive"));
00081   if (_length != length) {
00082     _length = length;
00083     NotifySubMeshesHypothesisModification();
00084   }
00085 }
00086 
00087 //=============================================================================
00091 //=============================================================================
00092 
00093 double StdMeshers_SegmentLengthAroundVertex::GetLength() const
00094 {
00095   return _length;
00096 }
00097 
00098 //=============================================================================
00102 //=============================================================================
00103 
00104 ostream & StdMeshers_SegmentLengthAroundVertex::SaveTo(ostream & save)
00105 {
00106   save << this->_length;
00107   return save;
00108 }
00109 
00110 //=============================================================================
00114 //=============================================================================
00115 
00116 istream & StdMeshers_SegmentLengthAroundVertex::LoadFrom(istream & load)
00117 {
00118   bool isOK = true;
00119   double a;
00120   isOK = (load >> a);
00121   if (isOK)
00122     this->_length = a;
00123   else
00124     load.clear(ios::badbit | load.rdstate());
00125   return load;
00126 }
00127 
00128 //=============================================================================
00132 //=============================================================================
00133 
00134 ostream & operator <<(ostream & save, StdMeshers_SegmentLengthAroundVertex & hyp)
00135 {
00136   return hyp.SaveTo( save );
00137 }
00138 
00139 //=============================================================================
00143 //=============================================================================
00144 
00145 istream & operator >>(istream & load, StdMeshers_SegmentLengthAroundVertex & hyp)
00146 {
00147   return hyp.LoadFrom( load );
00148 }
00149 
00150 //================================================================================
00157 //================================================================================
00158 
00159 bool StdMeshers_SegmentLengthAroundVertex::SetParametersByMesh(const SMESH_Mesh*   theMesh,
00160                                                                const TopoDS_Shape& theShape)
00161 {
00162   if ( !theMesh || theShape.IsNull() || theShape.ShapeType() != TopAbs_VERTEX )
00163     return false;
00164 
00165   SMESH_MeshEditor editor( const_cast<SMESH_Mesh*>( theMesh ) );
00166   SMESH_MesherHelper helper( *editor.GetMesh() );
00167 
00168   // get node built on theShape vertex
00169   SMESHDS_Mesh* meshDS = editor.GetMeshDS();
00170   SMESHDS_SubMesh* smV = meshDS->MeshElements( theShape );
00171   if ( !smV || smV->NbNodes() == 0 )
00172     return false;
00173   const SMDS_MeshNode* vNode = smV->GetNodes()->next();
00174 
00175   // calculate average length of segments sharing vNode
00176 
00177   _length = 0.;
00178   int nbSegs = 0;
00179 
00180   SMDS_ElemIteratorPtr segIt = vNode->GetInverseElementIterator(SMDSAbs_Edge);
00181   while ( segIt->more() ) {
00182     const SMDS_MeshElement* seg = segIt->next();
00183     // get geom edge
00184     int shapeID = editor.FindShape( seg );
00185     if (!shapeID) continue;
00186     const TopoDS_Shape& s = meshDS->IndexToShape( shapeID );
00187     if ( s.IsNull() || s.ShapeType() != TopAbs_EDGE ) continue;
00188     const TopoDS_Edge& edge = TopoDS::Edge( s );
00189     // params of edge ends
00190     double u0 = helper.GetNodeU( edge, seg->GetNode(0) );
00191     double u1 = helper.GetNodeU( edge, seg->GetNode(1) );
00192     // length
00193     BRepAdaptor_Curve AdaptCurve( edge );
00194     _length += GCPnts_AbscissaPoint::Length( AdaptCurve, u0, u1);
00195     nbSegs++;
00196   }
00197   
00198   if ( nbSegs > 1 )
00199     _length /= nbSegs;
00200 
00201   return nbSegs;
00202 }
00203 
00204 //================================================================================
00209 //================================================================================
00210 
00211 bool StdMeshers_SegmentLengthAroundVertex::SetParametersByDefaults(const TDefaults&,
00212                                                                    const SMESH_Mesh*)
00213 {
00214   return false;
00215 }
00216