From: ageay <ageay>
Date: Thu, 2 Feb 2012 10:28:00 +0000 (+0000)
Subject: *** empty log message ***
X-Git-Tag: V6_main_FINAL~890
X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=11730761a5e556b86e35fe3b80e13e4c5a0ef04f;p=tools%2Fmedcoupling.git

*** empty log message ***
---

diff --git a/src/MEDCoupling/MEDCouplingPointSet.cxx b/src/MEDCoupling/MEDCouplingPointSet.cxx
index a27ed5c3a..9c167c851 100644
--- a/src/MEDCoupling/MEDCouplingPointSet.cxx
+++ b/src/MEDCoupling/MEDCouplingPointSet.cxx
@@ -27,6 +27,7 @@
 #include "InterpKernelGeo2DQuadraticPolygon.hxx"
 #include "InterpKernelGeo2DNode.hxx"
 #include "DirectedBoundingBox.hxx"
+#include "InterpKernelAutoPtr.hxx"
 
 #include <cmath>
 #include <limits>
@@ -466,7 +467,7 @@ void MEDCouplingPointSet::tryToShareSameCoords(const MEDCouplingPointSet& other,
 void MEDCouplingPointSet::findNodesOnPlane(const double *pt, const double *vec, double eps, std::vector<int>& nodes) const throw(INTERP_KERNEL::Exception)
 {
   if(getSpaceDimension()!=3)
-    throw INTERP_KERNEL::Exception("Invalid spacedim to be applied on this ! Must be equal to 3 !");
+    throw INTERP_KERNEL::Exception("MEDCouplingPointSet::findNodesOnPlane : Invalid spacedim to be applied on this ! Must be equal to 3 !");
   int nbOfNodes=getNumberOfNodes();
   double a=vec[0],b=vec[1],c=vec[2],d=-pt[0]*vec[0]-pt[1]*vec[1]-pt[2]*vec[2];
   double deno=sqrt(a*a+b*b+c*c);
@@ -479,6 +480,55 @@ void MEDCouplingPointSet::findNodesOnPlane(const double *pt, const double *vec,
     }
 }
 
+/*!
+ * This method is expecting to be called for meshes so that getSpaceDimension() returns 2 or 3.
+ * This method returns in 'nodes' output all the nodes that are at a distance lower than epsilon from a line
+ * defined by the point 'pt' and the vector 'vec'. 'pt' and 'vec' are expected to have a dimension equal to space dimension of 'this'
+ * @param pt points to an array of size this->getSpaceDimension and represents a point that owns to plane.
+ * @param vec points to an array of size this->getSpaceDimension and represents the direction vector of the line. The norm of the vector is not compulsory equal to 1.
+ *            But norm must be greater than 10*abs(eps)
+ * @param eps is the maximal distance around the plane where node in this->_coords will be picked.
+ * @param nodes is the output of the method. The vector is not compulsory empty before call. The nodes that fulfills the condition will be added at the end of the nodes.
+ */
+void MEDCouplingPointSet::findNodesOnLine(const double *pt, const double *vec, double eps, std::vector<int>& nodes) const throw(INTERP_KERNEL::Exception)
+{
+  int spaceDim=getSpaceDimension();
+  if(spaceDim!=2 && spaceDim!=3)
+    throw INTERP_KERNEL::Exception("MEDCouplingPointSet::findNodesOnLine : Invalid spacedim to be applied on this ! Must be equal to 2 or 3 !");
+  int nbOfNodes=getNumberOfNodes();
+  int den=0.;
+  for(int i=0;i<spaceDim;i++)
+    den+=vec[i]*vec[i];
+  double deno=sqrt(den);
+  if(deno<10.*eps)
+    throw INTERP_KERNEL::Exception("MEDCouplingPointSet::findNodesOnLine : Invalid given direction vector ! Norm is too small !");
+  INTERP_KERNEL::AutoPtr<double> vecn=new double[spaceDim];
+  for(int i=0;i<spaceDim;i++)
+    vecn[i]=vec[i]/deno;
+  const double *work=_coords->getConstPointer();
+  if(spaceDim==2)
+    {
+      for(int i=0;i<nbOfNodes;i++)
+        {
+          if(std::abs(vecn[0]*(work[1]-pt[1])-vecn[1]*(work[0]-pt[0]))<eps)
+            nodes.push_back(i);
+          work+=2;
+        }
+    }
+  else
+    {
+      for(int i=0;i<nbOfNodes;i++)
+        {
+          double a=vecn[0]*(work[1]-pt[1])-vecn[1]*(work[0]-pt[0]);
+          double b=vecn[1]*(work[2]-pt[2])-vecn[2]*(work[1]-pt[1]);
+          double c=vecn[2]*(work[0]-pt[0])-vecn[0]*(work[2]-pt[2]);
+          if(std::sqrt(a*a+b*b+c*c)<eps)
+            nodes.push_back(i);
+          work+=3;
+        }
+    }
+}
+
 /*!
  * merge _coords arrays of m1 and m2 and returns the union. The returned instance is newly created with ref count == 1.
  */
diff --git a/src/MEDCoupling/MEDCouplingPointSet.hxx b/src/MEDCoupling/MEDCouplingPointSet.hxx
index d446775e3..e84d46fa5 100644
--- a/src/MEDCoupling/MEDCouplingPointSet.hxx
+++ b/src/MEDCoupling/MEDCouplingPointSet.hxx
@@ -80,6 +80,7 @@ namespace ParaMEDMEM
     void tryToShareSameCoords(const MEDCouplingPointSet& other, double epsilon) throw(INTERP_KERNEL::Exception);
     virtual void tryToShareSameCoordsPermute(const MEDCouplingPointSet& other, double epsilon) throw(INTERP_KERNEL::Exception) = 0;
     void findNodesOnPlane(const double *pt, const double *vec, double eps, std::vector<int>& nodes) const throw(INTERP_KERNEL::Exception);
+    void findNodesOnLine(const double *pt, const double *vec, double eps, std::vector<int>& nodes) const throw(INTERP_KERNEL::Exception);
     static DataArrayDouble *MergeNodesArray(const MEDCouplingPointSet *m1, const MEDCouplingPointSet *m2) throw(INTERP_KERNEL::Exception);
     static DataArrayDouble *MergeNodesArray(const std::vector<const MEDCouplingPointSet *>& ms) throw(INTERP_KERNEL::Exception);
     static MEDCouplingPointSet *BuildInstanceFromMeshType(MEDCouplingMeshType type);
diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx
index 2eb8fee7e..3f7fb4ffb 100644
--- a/src/MEDCoupling/MEDCouplingUMesh.cxx
+++ b/src/MEDCoupling/MEDCouplingUMesh.cxx
@@ -4843,6 +4843,8 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo
       std::map<int,INTERP_KERNEL::Node *> mappRev;
       INTERP_KERNEL::QuadraticPolygon pol1;
       MEDCouplingUMeshBuildQPFromMesh3(coo1,offset1,coo2,offset2,addCoords,b1[i],desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,/* output */pol1,mapp,mappRev);
+      std::vector<int> crTmp,crITmp;
+      crITmp.push_back(crI.back());
       for(std::vector<int>::const_iterator it2=candidates2.begin();it2!=candidates2.end();it2++)
         {
           INTERP_KERNEL::QuadraticPolygon pol2;
@@ -4850,6 +4852,11 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo
           MEDCouplingUMeshBuildQPFromMesh3(coo1,offset1,coo2,offset2,addCoords,b2[*it2],desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2,/* output */pol2,mapp,mappRev);
           pol1.buildPartitionsAbs(pol2,mapp,i,*it2,cr,crI,cNb1,cNb2);
         }
+      if(!crTmp.empty())
+        {
+          cr.insert(cr.end(),crTmp.begin(),crTmp.end());
+          crI.insert(crI.end(),crITmp.begin()+1,crITmp.end());
+        }
       for(std::map<int,INTERP_KERNEL::Node *>::const_iterator it=mappRev.begin();it!=mappRev.end();it++)
         (*it).second->decrRef();
     }
@@ -4961,6 +4968,13 @@ void MEDCouplingUMesh::BuildIntersectEdges(const MEDCouplingUMesh *m1, const MED
     }
 }
 
+/*!
+ * Given a 2D mesh conn by (conn2D,connI2D) it returns a single polygon
+ */
+void MEDCouplingUMesh::BuildUnionOf2DMesh(const std::vector<int>& conn2D, const std::vector<int>& connI2D, std::vector<int>& polyUnion)
+{
+}
+
 MEDCouplingUMeshCellIterator::MEDCouplingUMeshCellIterator(MEDCouplingUMesh *mesh):_mesh(mesh),_cell(new MEDCouplingUMeshCell(mesh)),
                                                                                    _own_cell(true),_cell_id(-1),_nb_cell(0)
 {
diff --git a/src/MEDCoupling/MEDCouplingUMesh.hxx b/src/MEDCoupling/MEDCouplingUMesh.hxx
index e5ed8196a..5f2bbadb1 100644
--- a/src/MEDCoupling/MEDCouplingUMesh.hxx
+++ b/src/MEDCoupling/MEDCouplingUMesh.hxx
@@ -225,6 +225,7 @@ namespace ParaMEDMEM
                                                   const MEDCouplingUMesh *m2, const std::vector<bool>& b2, const int *desc2, const int *descIndx2, const std::vector<std::vector<int> >& intesctEdges2,
                                                   const std::vector<double>& addCoords,
                                                   std::vector<int>& cr, std::vector<int>& crI, std::vector<int>& cNb1, std::vector<int>& cNb2);
+    static void BuildUnionOf2DMesh(const std::vector<int>& conn2D, const std::vector<int>& connI2D, std::vector<int>& polyUnion);
 /// @endcond
   private:
     //! this iterator stores current position in _nodal_connec array.
diff --git a/src/MEDCoupling/Test/MEDCouplingBasicsTest4.cxx b/src/MEDCoupling/Test/MEDCouplingBasicsTest4.cxx
index 5d7fc310c..bbc255b7b 100644
--- a/src/MEDCoupling/Test/MEDCouplingBasicsTest4.cxx
+++ b/src/MEDCoupling/Test/MEDCouplingBasicsTest4.cxx
@@ -2036,3 +2036,28 @@ void MEDCouplingBasicsTest4::testIntersect2DMeshesTmp1()
   coordY->decrRef();
   m1c->decrRef();
 }
+
+void MEDCouplingBasicsTest4::testFindNodesOnLine1()
+{
+  MEDCouplingUMesh *mesh=build2DTargetMesh_1();
+  const double pt[2]={-0.3,-0.3};
+  const double pt2[3]={0.,0.,0.};
+  const double pt3[3]={-0.3,0.,0.};
+  const double vec[2]={0.,1.};
+  const double vec2[3]={1.,0.,0.};
+  const double vec3[3]={0.,1.,1.};
+  const int expected1[3]={0,3,6};
+  std::vector<int> res;
+  mesh->findNodesOnLine(pt,vec,1e-12,res);
+  CPPUNIT_ASSERT_EQUAL(3,(int)res.size());
+  CPPUNIT_ASSERT(std::equal(expected1,expected1+3,res.begin()));
+  res.clear();
+  //
+  mesh->changeSpaceDimension(3);
+  mesh->rotate(pt2,vec2,M_PI/4.);
+  mesh->findNodesOnLine(pt3,vec3,1e-12,res);
+  CPPUNIT_ASSERT_EQUAL(3,(int)res.size());
+  CPPUNIT_ASSERT(std::equal(expected1,expected1+3,res.begin()));
+  //
+  mesh->decrRef();
+}
diff --git a/src/MEDCoupling/Test/MEDCouplingBasicsTest4.hxx b/src/MEDCoupling/Test/MEDCouplingBasicsTest4.hxx
index bcb45a954..d36129c98 100644
--- a/src/MEDCoupling/Test/MEDCouplingBasicsTest4.hxx
+++ b/src/MEDCoupling/Test/MEDCouplingBasicsTest4.hxx
@@ -89,6 +89,7 @@ namespace ParaMEDMEM
     CPPUNIT_TEST( testGetNodeIdsInUse1 );
     CPPUNIT_TEST( testBuildDescendingConnec2 );
     CPPUNIT_TEST( testIntersect2DMeshesTmp1 );
+    CPPUNIT_TEST( testFindNodesOnLine1 );
     CPPUNIT_TEST_SUITE_END();
   public:
     void testDescriptionInMeshTimeUnit1();
@@ -147,6 +148,7 @@ namespace ParaMEDMEM
     void testGetNodeIdsInUse1();
     void testBuildDescendingConnec2();
     void testIntersect2DMeshesTmp1();
+    void testFindNodesOnLine1();
   };
 }
 
diff --git a/src/MEDCoupling_Swig/MEDCoupling.i b/src/MEDCoupling_Swig/MEDCoupling.i
index 5ccd61b35..60b2f3d5c 100644
--- a/src/MEDCoupling_Swig/MEDCoupling.i
+++ b/src/MEDCoupling_Swig/MEDCoupling.i
@@ -817,6 +817,20 @@ namespace ParaMEDMEM
                  self->renumberNodes2(da2->getConstPointer(),newNbOfNodes);
                }
            }
+           PyObject *findNodesOnLine(PyObject *pt, PyObject *vec, double eps) const throw(INTERP_KERNEL::Exception)
+             {
+               std::vector<int> nodes;
+               int spaceDim=self->getSpaceDimension();
+               int sz1,sz2;
+               double *p=convertPyToNewDblArr2(pt,&sz1);
+               double *v=convertPyToNewDblArr2(vec,&sz2);
+               if(sz1!=spaceDim || sz2!=spaceDim)
+                 throw INTERP_KERNEL::Exception("Mismatch of spaceDimension and the length of the input array point and vector !");
+               self->findNodesOnLine(p,v,eps,nodes);
+               delete [] v;
+               delete [] p;
+               return convertIntArrToPyList2(nodes);
+             }
            PyObject *findNodesOnPlane(PyObject *pt, PyObject *vec, double eps) const throw(INTERP_KERNEL::Exception)
              {
                std::vector<int> nodes;
diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py
index 36e41b858..0091508ea 100644
--- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py
+++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest.py
@@ -8130,6 +8130,26 @@ class MEDCouplingBasicsTest(unittest.TestCase):
             self.assertAlmostEqual(expected5[i],m3.getCoords().getIJ(0,i),12);
             pass
         pass
+
+    def testFindNodesOnLine1(self):
+        mesh=MEDCouplingDataForTest.build2DTargetMesh_1();
+        pt=[-0.3,-0.3]
+        pt2=[0.,0.,0.]
+        pt3=[-0.3,0.,0.]
+        vec=[0.,1.]
+        vec2=[1.,0.,0.]
+        vec3=[0.,1.,1.]
+        expected1=[0,3,6]
+        res=mesh.findNodesOnLine(pt,vec,1e-12);
+        self.assertEqual(3,len(res));
+        self.assertEqual(expected1,res);
+        #
+        mesh.changeSpaceDimension(3);
+        mesh.rotate(pt2,vec2,pi/4.);
+        res=mesh.findNodesOnLine(pt3,vec3,1e-12);
+        self.assertEqual(3,len(res));
+        self.assertEqual(expected1,res);
+        pass
     
     def setUp(self):
         pass