X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FMEDCoupling%2FMEDCouplingMemArray.cxx;h=53213b741d529f377d1525d69a924f4e14631e14;hb=b7e277ad50223814bc479c5bd64cf12abe8fc959;hp=c70eba32ec646ee8250d64ad6ffb37e3f58b3041;hpb=93e2174a92103f402b8b68a75a3124c89d16cd7e;p=tools%2Fmedcoupling.git

diff --git a/src/MEDCoupling/MEDCouplingMemArray.cxx b/src/MEDCoupling/MEDCouplingMemArray.cxx
index c70eba32e..53213b741 100644
--- a/src/MEDCoupling/MEDCouplingMemArray.cxx
+++ b/src/MEDCoupling/MEDCouplingMemArray.cxx
@@ -36,6 +36,8 @@ typedef double (*MYFUNCPTR)(double);
 
 using namespace MEDCoupling;
 
+template class MemArray<int>;
+template class MemArray<double>;
 template class DataArrayTemplate<int>;
 template class DataArrayTemplate<double>;
 
@@ -112,6 +114,22 @@ void DataArrayDouble::FindClosestTupleIdAlg(const BBTreePts<SPACEDIM,int>& myTre
     }
 }
 
+int DataArray::EffectiveCircPerm(int nbOfShift, int nbOfTuples)
+{
+  if(nbOfTuples<=0)
+    throw INTERP_KERNEL::Exception("DataArray::EffectiveCircPerm : number of tuples is expected to be > 0 !");
+  if(nbOfShift>=0)
+    {
+      return nbOfShift%nbOfTuples;
+    }
+  else
+    {
+      int tmp(-nbOfShift);
+      tmp=tmp%nbOfTuples;
+      return nbOfTuples-tmp;
+    }
+}
+
 std::size_t DataArray::getHeapMemorySizeWithoutChildren() const
 {
   std::size_t sz1=_name.capacity();
@@ -2046,6 +2064,7 @@ DataArrayDouble *DataArrayDouble::accumulatePerChunck(const int *bgOfIndex, cons
  *          is to delete this array using decrRef() as it is no more needed. The array
  *          does not contain any textual info on components.
  *  \throw If \a this->getNumberOfComponents() != 2.
+ * \sa fromCartToPolar
  */
 DataArrayDouble *DataArrayDouble::fromPolarToCart() const
 {
@@ -2076,6 +2095,7 @@ DataArrayDouble *DataArrayDouble::fromPolarToCart() const
  *          on the third component is copied from \a this array. The caller
  *          is to delete this array using decrRef() as it is no more needed. 
  *  \throw If \a this->getNumberOfComponents() != 3.
+ * \sa fromCartToCyl
  */
 DataArrayDouble *DataArrayDouble::fromCylToCart() const
 {
@@ -2108,6 +2128,7 @@ DataArrayDouble *DataArrayDouble::fromCylToCart() const
  *          on the third component is copied from \a this array. The caller
  *          is to delete this array using decrRef() as it is no more needed.
  *  \throw If \a this->getNumberOfComponents() != 3.
+ * \sa fromCartToSpher
  */
 DataArrayDouble *DataArrayDouble::fromSpherToCart() const
 {
@@ -2179,6 +2200,116 @@ DataArrayDouble *DataArrayDouble::cartesianize(MEDCouplingAxisType atOfThis) con
   return ret.retn();
 }
 
+/*!
+ * This method returns a newly created array to be deallocated that contains the result of conversion from cartesian to polar.
+ * This method expects that \a this has exactly 2 components.
+ * \sa fromPolarToCart
+ */
+DataArrayDouble *DataArrayDouble::fromCartToPolar() const
+{
+  MCAuto<DataArrayDouble> ret(DataArrayDouble::New());
+  checkAllocated();
+  int nbOfComp(getNumberOfComponents()),nbTuples(getNumberOfTuples());
+  if(nbOfComp!=2)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::fromCartToPolar : must be an array with exactly 2 components !");
+  ret->alloc(nbTuples,2);
+  double *retPtr(ret->getPointer());
+  const double *ptr(begin());
+  for(int i=0;i<nbTuples;i++,ptr+=2,retPtr+=2)
+    {
+      retPtr[0]=sqrt(ptr[0]*ptr[0]+ptr[1]*ptr[1]);
+      retPtr[1]=atan2(ptr[1],ptr[0]);
+    }
+  return ret.retn();
+}
+
+/*!
+ * This method returns a newly created array to be deallocated that contains the result of conversion from cartesian to cylindrical.
+ * This method expects that \a this has exactly 3 components.
+ * \sa fromCylToCart
+ */
+DataArrayDouble *DataArrayDouble::fromCartToCyl() const
+{
+  MCAuto<DataArrayDouble> ret(DataArrayDouble::New());
+  checkAllocated();
+  int nbOfComp(getNumberOfComponents()),nbTuples(getNumberOfTuples());
+  if(nbOfComp!=3)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::fromCartToCyl : must be an array with exactly 3 components !");
+  ret->alloc(nbTuples,3);
+  double *retPtr(ret->getPointer());
+  const double *ptr(begin());
+  for(int i=0;i<nbTuples;i++,ptr+=3,retPtr+=3)
+    {
+      retPtr[0]=sqrt(ptr[0]*ptr[0]+ptr[1]*ptr[1]);
+      retPtr[1]=atan2(ptr[1],ptr[0]);
+      retPtr[2]=ptr[2];
+    }
+  return ret.retn();
+}
+
+/*!
+ * This method returns a newly created array to be deallocated that contains the result of conversion from cartesian to spherical coordinates.
+ * \sa fromSpherToCart
+ */
+DataArrayDouble *DataArrayDouble::fromCartToSpher() const
+{
+  MCAuto<DataArrayDouble> ret(DataArrayDouble::New());
+  checkAllocated();
+  int nbOfComp(getNumberOfComponents()),nbTuples(getNumberOfTuples());
+  if(nbOfComp!=3)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::fromCartToSpher : must be an array with exactly 3 components !");
+  ret->alloc(nbTuples,3);
+  double *retPtr(ret->getPointer());
+  const double *ptr(begin());
+  for(int i=0;i<nbTuples;i++,ptr+=3,retPtr+=3)
+    {
+      retPtr[0]=sqrt(ptr[0]*ptr[0]+ptr[1]*ptr[1]+ptr[2]*ptr[2]);
+      retPtr[1]=acos(ptr[2]/retPtr[0]);
+      retPtr[2]=atan2(ptr[1],ptr[0]);
+    }
+  return ret.retn();
+}
+
+/*!
+ * This method returns a newly created array to be deallocated that contains the result of conversion from cartesian to cylindrical relative to the given \a center and a \a vector.
+ * This method expects that \a this has exactly 3 components.
+ * \sa MEDCouplingFieldDouble::computeVectorFieldCyl
+ */
+DataArrayDouble *DataArrayDouble::fromCartToCylGiven(const DataArrayDouble *coords, const double center[3], const double vect[3]) const
+{
+  if(!coords)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::fromCartToCylGiven : input coords are NULL !");
+  MCAuto<DataArrayDouble> ret(DataArrayDouble::New());
+  checkAllocated(); coords->checkAllocated();
+  int nbOfComp(getNumberOfComponents()),nbTuples(getNumberOfTuples());
+  if(nbOfComp!=3)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::fromCartToCylGiven : must be an array with exactly 3 components !");
+  if(coords->getNumberOfComponents()!=3)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::fromCartToCylGiven : coords array must have exactly 3 components !");
+  if(coords->getNumberOfTuples()!=nbTuples)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::fromCartToCylGiven : coords array must have the same number of tuples !");
+  ret->alloc(nbTuples,nbOfComp);
+  double magOfVect(sqrt(vect[0]*vect[0]+vect[1]*vect[1]+vect[2]*vect[2]));
+  if(magOfVect<1e-12)
+    throw INTERP_KERNEL::Exception("DataArrayDouble::fromCartToCylGiven : magnitude of vect is too low !");
+  double Ur[3],Uteta[3],Uz[3],*retPtr(ret->getPointer());
+  const double *coo(coords->begin()),*vectField(begin());
+  std::transform(vect,vect+3,Uz,std::bind2nd(std::multiplies<double>(),1./magOfVect));
+  for(int i=0;i<nbTuples;i++,vectField+=3,retPtr+=3,coo+=3)
+    {
+      std::transform(coo,coo+3,center,Ur,std::minus<double>());
+      Uteta[0]=Uz[1]*Ur[2]-Uz[2]*Ur[1]; Uteta[1]=Uz[2]*Ur[0]-Uz[0]*Ur[2]; Uteta[2]=Uz[0]*Ur[1]-Uz[1]*Ur[0];
+      double magOfTeta(sqrt(Uteta[0]*Uteta[0]+Uteta[1]*Uteta[1]+Uteta[2]*Uteta[2]));
+      std::transform(Uteta,Uteta+3,Uteta,std::bind2nd(std::multiplies<double>(),1./magOfTeta));
+      Ur[0]=Uteta[1]*Uz[2]-Uteta[2]*Uz[1]; Ur[1]=Uteta[2]*Uz[0]-Uteta[0]*Uz[2]; Ur[2]=Uteta[0]*Uz[1]-Uteta[1]*Uz[0];
+      retPtr[0]=Ur[0]*vectField[0]+Ur[1]*vectField[1]+Ur[2]*vectField[2];
+      retPtr[1]=Uteta[0]*vectField[0]+Uteta[1]*vectField[1]+Uteta[2]*vectField[2];
+      retPtr[2]=Uz[0]*vectField[0]+Uz[1]*vectField[1]+Uz[2]*vectField[2];
+    }
+  ret->copyStringInfoFrom(*this);
+  return ret.retn();
+}
+
 /*!
  * Computes the doubly contracted product of every tensor defined by the tuple of \a this
  * array contating 6 components.
@@ -5507,6 +5638,48 @@ DataArrayInt *DataArrayInt::buildPermutationArr(const DataArrayInt& other) const
   return ret.retn();
 }
 
+/*!
+ * Elements of \a partOfThis are expected to be included in \a this.
+ * The returned array \a ret is so that this[ret]==partOfThis
+ *
+ * For example, if \a this array contents are [9,10,0,6,4,11,3,8] and if \a partOfThis contains [6,0,11,8]
+ * the return array will contain [3,2,5,7].
+ *
+ * \a this is expected to be a 1 compo allocated array.
+ * \param [in] partOfThis - A 1 compo allocated array
+ * \return - A newly allocated array to be dealed by caller having the same number of tuples than \a partOfThis.
+ * \throw if two same element is present twice in \a this
+ * \throw if an element in \a partOfThis is \b NOT in \a this.
+ */
+DataArrayInt *DataArrayInt::indicesOfSubPart(const DataArrayInt& partOfThis) const
+{
+  if(getNumberOfComponents()!=1 || partOfThis.getNumberOfComponents()!=1)
+    throw INTERP_KERNEL::Exception("DataArrayInt::indicesOfSubPart : this and input array must be one component array !");
+  checkAllocated(); partOfThis.checkAllocated();
+  int thisNbTuples(getNumberOfTuples()),nbTuples(partOfThis.getNumberOfTuples());
+  const int *thisPt(begin()),*pt(partOfThis.begin());
+  MCAuto<DataArrayInt> ret(DataArrayInt::New());
+  ret->alloc(nbTuples,1);
+  int *retPt(ret->getPointer());
+  std::map<int,int> m;
+  for(int i=0;i<thisNbTuples;i++,thisPt++)
+    m[*thisPt]=i;
+  if(m.size()!=thisNbTuples)
+    throw INTERP_KERNEL::Exception("DataArrayInt::indicesOfSubPart : some elements appears more than once !");
+  for(int i=0;i<nbTuples;i++,retPt++,pt++)
+    {
+      std::map<int,int>::const_iterator it(m.find(*pt));
+      if(it!=m.end())
+        *retPt=(*it).second;
+      else
+        {
+          std::ostringstream oss; oss << "DataArrayInt::indicesOfSubPart : At pos #" << i << " of input array value is " << *pt << " not in this !";
+          throw INTERP_KERNEL::Exception(oss.str());
+        }
+    }
+  return ret.retn();
+}
+
 void DataArrayInt::aggregate(const DataArrayInt *other)
 {
   if(!other)