Changeset 1390 in Sophya

Timestamp:

Feb 9, 2001, 6:08:57 PM (25 years ago)

Author:

ansari

Message:

Modifs interface FFTServer , Reza 9/2/2001

Location:

trunk/SophyaLib/NTools

Files:

• cimage.cc (modified) (5 diffs)
• cimage.h (modified) (7 diffs)
• fftpserver.cc (modified) (3 diffs)
• fftpserver.h (modified) (4 diffs)
• fftservintf.cc (modified) (3 diffs)
• fftservintf.h (modified) (2 diffs)

trunk/SophyaLib/NTools/cimage.cc

@@ -9 +9 @@
 #include "cimage.h"
 
+/*!
+   \class SOPHYA::Image
+   \ingroup NTools
+   This class specializes the Matrix class, for representing intensity 
+   (or grey-level) images. It adds the possibility of defining pixel
+   size and offset. The convention for the acces operator is also
+   changed, compared to the matrix class.
+*/
+
+//! Default constructor - Pixel array not allocated
 template <class T>
 Image<T>::Image()
@@ -18 +28 @@
 }
 
-template <class T>
-Image<T>::Image(uint_4 sizx, uint_4 sizy, r_8 szpx, r_8 szpy, r_8 orgx, r_8 orgy)
+//! Constructor, with specification of the image size, pixel size and offsets
+template <class T>
+Image<T>::Image(sa_size_t sizx, sa_size_t sizy, r_8 szpx, r_8 szpy, r_8 orgx, r_8 orgy)
 // Constructeur par defaut.
   : TMatrix<T>(sizy, sizx) 
@@ -27 +38 @@
 }
 
+//! Copy constructor
 template <class T>
 Image<T>::Image(const Image<T>& a)
@@ -35 +47 @@
 }
 
+
+//! Copy constructor with possibility of duplicating the pixel array
+/*!
+  \param share : if true, the data is shared, duplicated if false.
+*/
 template <class T>
 Image<T>::Image(const Image<T>& a, bool share)
@@ -42 +59 @@
   SetPixelSize(a.XPixSize(), a.YPixSize());
 }
+
+//! Destructor 
 template <class T>
 Image<T>::~Image()

trunk/SophyaLib/NTools/cimage.h

@@ -16 +16 @@
 namespace SOPHYA {
 
+//! Class for handling images
 template <class T>
 class Image : public TMatrix<T> {
@@ -21 +22 @@
 public:
   Image();
-  Image(uint_4 sizx, uint_4 sizy, r_8 szpx=1., r_8 szpy=1., r_8 orgx=0., r_8 orgy=0.);
+  Image(sa_size_t sizx, sa_size_t sizy, r_8 szpx=1., r_8 szpy=1., r_8 orgx=0., r_8 orgy=0.);
   Image(const Image<T>& a);
   Image(const Image<T>& a, bool share);
@@ -28 +29 @@
 
   // Inline element acces methods 
-  inline T const& operator()(uint_4 ix,uint_4 jy) const;
-  inline T&       operator()(uint_4 ix,uint_4 jy);
+  inline T const& operator()(sa_size_t ix, sa_size_t jy) const;
+  inline T&       operator()(sa_size_t ix, sa_size_t jy);
 
-  inline uint_4  XSize() const {return NCols();}
-  inline uint_4  YSize() const {return NRows();}
+  //! Returns the image size along X (corresponding to the number of columns)  
+  inline sa_size_t  XSize() const {return NCols();}
+  //! Returns the image size along Y (corresponding to the number of lines)  
+  inline sa_size_t  YSize() const {return NRows();}
 
+  //! Returns the X position, for pixel(0,0)
   inline r_8     XOrg() const {return org_x;}
+  //! Returns the Y position, for pixel(0,0)
   inline r_8     YOrg() const {return org_y;}
 
+  //! Returns pixel size along X 
   inline r_8     XPixSize() const {return pxsz_x;}
+  //! Returns pixel size along Y
   inline r_8     YPixSize() const {return pxsz_y;}
 
-  inline r_8     XPos(uint_4 i) const {return org_x + i*pxsz_x;}
-  inline r_8     YPos(uint_4 j) const {return org_y + j*pxsz_y;}
+  //! Returns the pixel position along X, for pixels with X index \b i
+  inline r_8     XPos(sa_size_t i) const {return org_x + i*pxsz_x;}
+  //! Returns the pixel position along Y, for pixels with Y index \b j
+  inline r_8     YPos(sa_size_t j) const {return org_y + j*pxsz_y;}
 
+  //! Set the position for pixel(0,0)
   inline void    SetOrg(r_8 orgx=0., r_8 orgy=0.)
                        { org_x = orgx;  org_y = orgy; }  
+  //! Set the pixel size
   inline void    SetPixelSize(r_8 szx=1., r_8 szy=1.)
                        { pxsz_x = szx;   pxsz_y = szy; }
@@ -53 +64 @@
 };
 
+//! () : Return the pixel value (element) for pixel \b ix (column index) and \b jy (line index)
 template <class T>
-inline T const& Image<T>::operator()(uint_4 ix, uint_4 jy) const
+inline T const& Image<T>::operator()(sa_size_t ix, sa_size_t jy) const
 {
   return (TMatrix<T>::operator() (jy, ix));
 }
 
+//! () : Return the pixel value (element) for pixel \b ix (column index) and \b jy (line index)
 template <class T>
-inline T & Image<T>::operator()(uint_4 ix, uint_4 jy) 
+inline T & Image<T>::operator()(sa_size_t ix, sa_size_t jy) 
 {
   return (TMatrix<T>::operator() (jy, ix));
@@ -82 +95 @@
 };
 
-/*! \ingroup Image \fn operator<<(POutPersist&,Image<T>&)
+/*! \ingroup NTools \fn operator<<(POutPersist&,Image<T>&)
   \brief Write Image \b obj into POutPersist stream \b os */
 template <class T>
@@ -88 +101 @@
 { FIO_Image<T> fio(&obj);  fio.Write(os);  return(os); }
 
-/*! \ingroup Image \fn operator>>(PInPersist&,Image<T>&)
+/*! \ingroup NTools \fn operator>>(PInPersist&,Image<T>&)
   \brief Read Image \b obj from PInPersist stream \b os */
 template <class T>
@@ -102 +115 @@
 
 class GeneralFit;
+//! OBSOLETE - Should not be used - for compatibility with peida in (s)piapp
 class RzImage {
 public:

trunk/SophyaLib/NTools/fftpserver.cc

@@ -28 +28 @@
 FFTPackServer::FFTPackServer()
    : FFTServerInterface("FFTPackServer using extended FFTPack (C-version) package")
-
+  , ckR4(true, true) , ckR8(true, true)
 {
 //the working array and its size for the different
@@ -56 +56 @@
 
 
-void FFTPackServer::FFTForward(TVector< complex<r_8> > const & in, TVector< complex<r_8> > & out)
-{
-  out = in;
-  fftf(out.NElts(), out.Data());
-  if (getNormalize()) out *= (1./(r_8)(in.NElts()));
-}
-
-void FFTPackServer::FFTBackward(TVector< complex<r_8> > const & in, TVector< complex<r_8> > & out)
-{
-  out = in;
-  fftb(out.NElts(), out.Data());
-}
-
-
-
-void FFTPackServer::FFTForward(TVector< complex<r_4> > const & in, TVector< complex<r_4> > & out)
-{
-  out = in;
-  fftf(out.NElts(), out.Data());
-  if (getNormalize()) out *= (1./(r_4)(in.NElts()));
-}
-
-void FFTPackServer::FFTBackward(TVector< complex<r_4> > const & in, TVector< complex<r_4> > & out)
-{
-  out = in;
-  fftb(out.NElts(), out.Data());
-}
-
-void FFTPackServer::FFTForward(TVector< r_4 > const & in, TVector< complex<r_4> > & out)
-{
-  TVector< r_4 > inout(in, false);
-  fftf(inout.NElts(), inout.Data());
-  ReShapetoCompl(inout, out);
-  if (getNormalize()) out *= (1./(r_4)(in.NElts()));
-}
-
-void FFTPackServer::FFTBackward(TVector< complex<r_4> > const & in, TVector< r_4 > & out)
-{
-  ReShapetoReal(in, out);
-  fftb(out.NElts(), out.Data());
-}
-
-
-void FFTPackServer::FFTForward(TVector< r_8 > const & in, TVector< complex<r_8> > & out)
-{
-  TVector< r_8 > inout(in, false);
-  fftf(inout.NElts(), inout.Data());
-  ReShapetoCompl(inout, out);
-  if (getNormalize()) out *= (1./(r_8)(in.NElts()));
-}
-
-void FFTPackServer::FFTBackward(TVector< complex<r_8> > const & in, TVector< r_8 > & out)
-{
-  ReShapetoReal(in, out);
-  fftb(out.NElts(), out.Data());
+void FFTPackServer::FFTForward(TArray< complex<r_8> > const & in, TArray< complex<r_8> > & out)
+{
+  ckR8.CheckResize(in, out);
+  out = in;
+  fftf(out.Size(), out.Data());
+  if (getNormalize()) out *= (1./(r_8)(in.Size()));
+}
+
+void FFTPackServer::FFTBackward(TArray< complex<r_8> > const & in, TArray< complex<r_8> > & out)
+{
+  ckR8.CheckResize(in, out);
+  out = in;
+  fftb(out.Size(), out.Data());
+}
+
+
+
+void FFTPackServer::FFTForward(TArray< complex<r_4> > const & in, TArray< complex<r_4> > & out)
+{
+  ckR4.CheckResize(in, out);  
+  out = in;
+  fftf(out.Size(), out.Data());
+  if (getNormalize()) out *= (1./(r_4)(in.Size()));
+}
+
+void FFTPackServer::FFTBackward(TArray< complex<r_4> > const & in, TArray< complex<r_4> > & out)
+{
+  ckR4.CheckResize(in, out);  
+  out = in;
+  fftb(out.Size(), out.Data());
+}
+
+void FFTPackServer::FFTForward(TArray< r_4 > const & in, TArray< complex<r_4> > & out)
+{
+  ckR4.CheckResize(in, out);  
+  TArray< r_4 > inout(in, false);
+  fftf(inout.Size(), inout.Data());
+  ckR4.ReShapetoCompl(inout, out);
+  if (getNormalize()) out *= complex<r_4>((1./(r_4)(in.Size())), 0.);
+}
+
+void FFTPackServer::FFTBackward(TArray< complex<r_4> > const & in, TArray< r_4 > & out)
+{
+  ckR4.CheckResize(in, out);    
+  ckR4.ReShapetoReal(in, out);
+  fftb(out.Size(), out.Data());
+}
+
+
+void FFTPackServer::FFTForward(TArray< r_8 > const & in, TArray< complex<r_8> > & out)
+{
+  ckR8.CheckResize(in, out);  
+  TArray< r_8 > inout(in, false);
+  fftf(inout.Size(), inout.Data());
+  ckR8.ReShapetoCompl(inout, out);
+  if (getNormalize()) out *= complex<r_8>((1./(r_8)(in.Size())), 0.);
+}
+
+void FFTPackServer::FFTBackward(TArray< complex<r_8> > const & in, TArray< r_8 > & out)
+{
+  ckR8.CheckResize(in, out);  
+  ckR8.ReShapetoReal(in, out);
+  fftb(out.Size(), out.Data());
 }
 
@@ -208 +216 @@
 }
 
-// Methodes pour reordonner les donnees 
-
-/* --Methode-- */
-void FFTPackServer::ReShapetoReal( TVector< complex<r_8> > const & in, TVector< r_8 >  & out)
-{
-  int n = in.NElts();
-  int ncs = (fabs(in(n-1).imag()) > 1.e-12) ? ncs = 2*n-1 : ncs = n*2-2; 
-  out.ReSize(ncs);
-  int k;
-  out(0) = in(0).real();
-  for(k=1;k<n-1;k++) {
-    out(2*k-1) = in(k).real();
-    out(2*k) = in(k).imag();
-  }
-  if (ncs == n*2-2)  out(ncs-1) = in(n-1).real();
-  else { out(ncs-2) = in(n-1).real();  out(ncs-1) = in(n-1).imag(); }
-}
-
-/* --Methode-- */
-void FFTPackServer::ReShapetoReal( TVector< complex<r_4> > const & in, TVector< r_4 >  & out)
-{
-  int n = in.NElts();
-  int ncs = (fabs(in(n-1).imag()) > 1.e-12) ? ncs = 2*n-1 : ncs = n*2-2; 
-  out.ReSize(ncs);
-  int k;
-  out(0) = in(0).real();
-  for(k=1;k<n-1;k++) {
-    out(2*k-1) = in(k).real();
-    out(2*k) = in(k).imag();
-  }
-  if (ncs == n*2-2)  out(ncs-1) = in(n-1).real();
-  else { out(ncs-2) = in(n-1).real();  out(ncs-1) = in(n-1).imag(); }
-}
-
-
-/* --Methode-- */
-void FFTPackServer::ReShapetoCompl(TVector< r_8 > const & in, TVector< complex<r_8> > & out)
-{
-  uint_4 n =  in.NElts();
-  uint_4 ncs = n/2+1;
-  uint_4 nc = (n%2 != 0) ? n/2+1 : n/2;
-  out.ReSize(ncs);
-  out(0) = complex<r_8> (in(0),0.);
-  for(int k=1;k<nc;k++) 
-    out(k) =  complex<r_4> (in(2*k-1), in(2*k));
-  if (n%2 == 0) out(ncs-1) = complex<r_8>(in(n-1), 0.);
-
-}
-
-/* --Methode-- */
-void FFTPackServer::ReShapetoCompl(TVector< r_4 > const & in, TVector< complex<r_4> > & out)
-{
-  uint_4 n =  in.NElts();
-  uint_4 ncs = n/2+1;
-  uint_4 nc = (n%2 != 0) ? n/2+1 : n/2;
-  out.ReSize(ncs);
-  out(0) = complex<r_4> (in(0),0.);
-  for(int k=1;k<nc;k++) 
-    out(k) =  complex<r_4> (in(2*k-1), in(2*k));
-  if (n%2 == 0) out(ncs-1) = complex<r_4>(in(n-1), 0.);
-}
+

trunk/SophyaLib/NTools/fftpserver.h

@@ -23 +23 @@
 
   // Transforme unidimensionnel sur des double
-  virtual void FFTForward(TVector< complex<r_8> > const & in, TVector< complex<r_8> > & out);
-  virtual void FFTBackward(TVector< complex<r_8> > const & in, TVector< complex<r_8> > & out);
-  virtual void FFTForward(TVector< r_8 > const & in, TVector< complex<r_8> > & out);
-  virtual void FFTBackward(TVector< complex<r_8> > const & in, TVector< r_8 > & out);
+  virtual void FFTForward(TArray< complex<r_8> > const & in, TArray< complex<r_8> > & out);
+  virtual void FFTBackward(TArray< complex<r_8> > const & in, TArray< complex<r_8> > & out);
+  virtual void FFTForward(TArray< r_8 > const & in, TArray< complex<r_8> > & out);
+  virtual void FFTBackward(TArray< complex<r_8> > const & in, TArray< r_8 > & out);
 
   // Transforme unidimensionnel sur des float
-  virtual void FFTForward(TVector< complex<r_4> > const & in, TVector< complex<r_4> > & out);
-  virtual void FFTBackward(TVector< complex<r_4> > const & in, TVector< complex<r_4> > & out);
-  virtual void FFTForward(TVector< r_4 > const & in, TVector< complex<r_4> > & out);
-  virtual void FFTBackward(TVector< complex<r_4> > const & in, TVector< r_4 > & out);
+  virtual void FFTForward(TArray< complex<r_4> > const & in, TArray< complex<r_4> > & out);
+  virtual void FFTBackward(TArray< complex<r_4> > const & in, TArray< complex<r_4> > & out);
+  virtual void FFTForward(TArray< r_4 > const & in, TArray< complex<r_4> > & out);
+  virtual void FFTBackward(TArray< complex<r_4> > const & in, TArray< r_4 > & out);
 
   //  Methodes propres a cette classe
@@ -43 +43 @@
   virtual void fftf(int_4 l, complex<double>* inout);
   virtual void fftb(int_4 l, complex<double>* inout);
-
-// Methodes statiques pour reordonner les donnees 
-  virtual  void ReShapetoReal( TVector< complex<r_8> > const & in, TVector< r_8 >  & out);
-  virtual  void ReShapetoReal( TVector< complex<r_4> > const & in, TVector< r_4 >  & out);
-
-  virtual void ReShapetoCompl(TVector< r_8 > const & in, TVector< complex<r_8> > & out);
-  virtual void ReShapetoCompl(TVector< r_4 > const & in, TVector< complex<r_4> > & out);
 
  protected:
@@ -68 +61 @@
   int sz_cdfft;
   r_8* ws_cdfft;
+
+  FFTArrayChecker<r_4> ckR4;
+  FFTArrayChecker<r_8> ckR8;
 };
 
@@ -73 +69 @@
 
 
-
-
 #endif

trunk/SophyaLib/NTools/fftservintf.cc

@@ -6 +6 @@
   \ingroup NTools
   Defines the interface for FFT (Fast Fourier Transform) operations.
-  Warning: The interface is likely to be changed and simplified,
-  using only SOPHYA::TArray Objects.
 */
 
@@ -22 +20 @@
 }
 
-/* --Methode-- */
-void FFTServerInterface::FFTForward(TVector< complex<r_8> > const &, TVector< complex<r_8> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTForward(TVector...) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTBackward(TVector< complex<r_8> > const &, TVector< complex<r_8> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTBackward(TVector...) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTForward(TVector< r_8 > const &, TVector< complex<r_8> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTForward(TVector...) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTBackward(TVector< complex<r_8> > const &, TVector< r_8 > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTBackward(TVector...) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTForward(TMatrix< complex<r_8> > const &, TMatrix< complex<r_8> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTForward(TMatrix...) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTBackward(TMatrix< complex<r_8> > const &, TMatrix< complex<r_8> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTBackward(TMatrix...) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTForward(TMatrix< r_8 > const &, TMatrix< complex<r_8> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTForward(TMatrix...) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTBackward(TMatrix< complex<r_8> > const &, TMatrix< r_8 > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTBackward(TMatrix...) Unsupported operation !");
+// ----------------- Transforme pour les double -------------------
+
+/* --Methode-- */
+void FFTServerInterface::FFTForward(TArray< complex<r_8> > const &, TArray< complex<r_8> > &)
+{
+  throw NotAvailableOperation("FFTServer::FFTForward(TArray...) Unsupported operation !");
+}
+
+/* --Methode-- */
+void FFTServerInterface::FFTBackward(TArray< complex<r_8> > const &, TArray< complex<r_8> > &)
+{
+  throw NotAvailableOperation("FFTServer::FFTBackward(TArray...) Unsupported operation !");
+}
+
+/* --Methode-- */
+void FFTServerInterface::FFTForward(TArray< r_8 > const &, TArray< complex<r_8> > &)
+{
+  throw NotAvailableOperation("FFTServer::FFTForward(TArray...) Unsupported operation !");
+}
+
+/* --Methode-- */
+void FFTServerInterface::FFTBackward(TArray< complex<r_8> > const &, TArray< r_8 > &)
+{
+  throw NotAvailableOperation("FFTServer::FFTBackward(TArray...) Unsupported operation !");
 }
 
@@ -74 +50 @@
 
 /* --Methode-- */
-void FFTServerInterface::FFTForward(TVector< complex<r_4> > const &, TVector< complex<r_4> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTForward(TVector r_4 ... ) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTBackward(TVector< complex<r_4> > const &, TVector< complex<r_4> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTBackward(TVector r_4 ... ) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTForward(TVector< r_4 > const &, TVector< complex<r_4> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTForward(TVector r_4 ... ) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTBackward(TVector< complex<r_4> > const &, TVector< r_4 > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTBackward(TVector r_4 ... ) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTForward(TMatrix< complex<r_4> > const &, TMatrix< complex<r_4> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTForward(TMatrix r_4 ... ) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTBackward(TMatrix< complex<r_4> > const &, TMatrix< complex<r_4> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTBackward(TMatrix r_4 ... ) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTForward(TMatrix< r_4 > const &, TMatrix< complex<r_4> > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTForward(TMatrix r_4 ... ) Unsupported operation !");
-}
-
-/* --Methode-- */
-void FFTServerInterface::FFTBackward(TMatrix< complex<r_4> > const &, TMatrix< r_4 > &)
-{
-  throw NotAvailableOperation("FFTServer::FFTBackward(TMatrix r_4 ... ) Unsupported operation !");
-}
-
-
+void FFTServerInterface::FFTForward(TArray< complex<r_4> > const &, TArray< complex<r_4> > &)
+{
+  throw NotAvailableOperation("FFTServer::FFTForward(TArray r_4 ... ) Unsupported operation !");
+}
+
+/* --Methode-- */
+void FFTServerInterface::FFTBackward(TArray< complex<r_4> > const &, TArray< complex<r_4> > &)
+{
+  throw NotAvailableOperation("FFTServer::FFTBackward(TArray r_4 ... ) Unsupported operation !");
+}
+
+/* --Methode-- */
+void FFTServerInterface::FFTForward(TArray< r_4 > const &, TArray< complex<r_4> > &)
+{
+  throw NotAvailableOperation("FFTServer::FFTForward(TArray r_4 ... ) Unsupported operation !");
+}
+
+/* --Methode-- */
+void FFTServerInterface::FFTBackward(TArray< complex<r_4> > const &, TArray< r_4 > &)
+{
+  throw NotAvailableOperation("FFTServer::FFTBackward(TArray r_4 ... ) Unsupported operation !");
+}
+
+
+
+/* --Methode-- */
+template <class T>
+FFTArrayChecker<T>::FFTArrayChecker(bool checkpack, bool onedonly)
+{
+  _checkpack = checkpack;
+  _onedonly = onedonly;
+}
+
+/* --Methode-- */
+template <class T>
+FFTArrayChecker<T>::~FFTArrayChecker()
+{
+}
+
+/* --Methode-- */
+template <class T>
+int FFTArrayChecker<T>::CheckResize(TArray< complex<T> > const & in, TArray< complex<T> > & out)
+{
+  int k;
+  if (in.Size() < 1) 
+    throw(SzMismatchError("FFTArrayChecker::CheckResize(complex in, complex out) - Unallocated input array !"));
+  if (_checkpack) 
+    if ( !in.IsPacked() )
+      throw(SzMismatchError("FFTArrayChecker::CheckResize(complex in, complex out) - Not packed input array !"));
+  int ndg1 = 0;
+  for(k=0; k<in.NbDimensions(); k++) 
+    if (in.Size(k) > 1)  ndg1++;
+  if (_onedonly) 
+    if (ndg1++ > 1)
+      throw(SzMismatchError("FFTArrayChecker::CheckResize(complex in, complex out) - Only 1-D array accepted !"));
+  
+  sa_size_t sz[BASEARRAY_MAXNDIMS];
+  for(k=0; k<in.NbDimensions(); k++) 
+    sz[k] = in.Size(k);
+  out.ReSize(in.NbDimensions(), sz);
+
+  return(ndg1);
+}
+
+/* --Methode-- */
+template <class T>
+int FFTArrayChecker<T>::CheckResize(TArray< T > const & in, TArray< complex<T> > & out)
+{
+  int k;
+  if (in.Size() < 1) 
+    throw(SzMismatchError("FFTArrayChecker::CheckResize(real in, complex out) - Unallocated input array !"));
+  if (_checkpack) 
+    if ( !in.IsPacked() )
+      throw(SzMismatchError("FFTArrayChecker::CheckResize(real in, complex out) - Not packed input array !"));
+  int ndg1 = 0;
+  for(k=0; k<in.NbDimensions(); k++) 
+    if (in.Size(k) > 1)  ndg1++;
+  if (_onedonly) 
+    if (ndg1++ > 1)
+      throw(SzMismatchError("FFTArrayChecker::CheckResize(real in, complex out) - Only 1-D array accepted !"));
+  
+  sa_size_t sz[BASEARRAY_MAXNDIMS];
+  for(k=0; k<in.NbDimensions(); k++) 
+    sz[k] = in.Size(k)/2+1;
+  out.ReSize(in.NbDimensions(), sz);
+
+  return(ndg1);
+}
+
+/* --Methode-- */
+template <class T>
+int FFTArrayChecker<T>::CheckResize(TArray< complex<T> > const & in, TArray< T > & out)
+{
+  int k;
+  if (in.Size() < 1) 
+    throw(SzMismatchError("FFTArrayChecker::CheckResize(complex in, real out) - Unallocated input array !"));
+  if (_checkpack) 
+    if ( !in.IsPacked() )
+      throw(SzMismatchError("FFTArrayChecker::CheckResize(complex in, real out) - Not packed input array !"));
+  int ndg1 = 0;
+  for(k=0; k<in.NbDimensions(); k++) 
+    if (in.Size(k) > 1)  ndg1++;
+  if (_onedonly) 
+    if (ndg1++ > 1)
+      throw(SzMismatchError("FFTArrayChecker::CheckResize(complex in, real out) - Only 1-D array accepted !"));
+  
+  sa_size_t sz[BASEARRAY_MAXNDIMS];
+  for(k=0; k<in.NbDimensions(); k++) 
+    sz[k] = in.Size(k)*2-1;
+  out.ReSize(in.NbDimensions(), sz);
+
+  return(ndg1);
+
+}
+
+/* --Methode-- */
+template <class T>
+void FFTArrayChecker<T>::ReShapetoReal( TArray< complex<T> > const & in, TArray< T >  & out)
+{
+  sa_size_t n = in.Size();
+  //  int ncs = (fabs(in(n-1).imag()) > 1.e-12) ? ncs = 2*n-1 : ncs = n*2-2;
+  sa_size_t ncs = 2*n-1;
+  sa_size_t k;
+  out[0] = in[0].real();
+  for(k=1;k<n-1;k++) {
+    out[2*k-1] = in[k].real();
+    out[2*k] = in[k].imag();
+  }
+  //  if (ncs == n*2-2)  out[ncs-1] = in[n-1].real();
+  //  else { out[ncs-2] = in[n-1].real();  out[ncs-1] = in[n-1].imag(); }
+  out[ncs-2] = in[n-1].real();  out[ncs-1] = in[n-1].imag(); 
+}
+
+/* --Methode-- */
+template <class T>
+void FFTArrayChecker<T>::ReShapetoCompl(TArray< T > const & in, TArray< complex<T> > & out)
+{
+  sa_size_t n =  in.Size();
+  sa_size_t ncs = n/2+1;
+  sa_size_t nc = (n%2 != 0) ? n/2+1 : n/2;
+  out[0] = complex<T> (in[0],0.);
+  for(int k=1;k<nc;k++) 
+    out[k] =  complex<T> (in[2*k-1], in[2*k]);
+  if (n%2 == 0) out[ncs-1] = complex<r_8>(in[n-1], 0.);
+}
+
+
+#ifdef __CXX_PRAGMA_TEMPLATES__
+#pragma define_template FFTArrayChecker<r_4>
+#pragma define_template FFTArrayChecker<r_8>
+#endif
+
+#if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES)
+template class FFTArrayChecker<r_4>;
+template class FFTArrayChecker<r_8>;
+#endif

trunk/SophyaLib/NTools/fftservintf.h

@@ -45 +45 @@
 
   //---------------------------------------------------
-  // Transforme 2D sur des doubles
-  virtual void FFTForward(TMatrix< complex<r_8> > const & in, TMatrix< complex<r_8> > & out);
-  virtual void FFTBackward(TMatrix< complex<r_8> > const & in, TMatrix< complex<r_8> > & out);
-  virtual void FFTForward(TMatrix< r_8 > const & in, TMatrix< complex<r_8> > & out);
-  virtual void FFTBackward(TMatrix< complex<r_8> > const & in, TMatrix< r_8 > & out);
+  // Transforme N-dim sur des doubles
+  virtual void FFTForward(TArray< complex<r_8> > const & in, TArray< complex<r_8> > & out);
+  virtual void FFTBackward(TArray< complex<r_8> > const & in, TArray< complex<r_8> > & out);
+  virtual void FFTForward(TArray< r_8 > const & in, TArray< complex<r_8> > & out);
+  virtual void FFTBackward(TArray< complex<r_8> > const & in, TArray< r_8 > & out);
 
-  // Transforme 2D sur des float
-  virtual void FFTForward(TMatrix< complex<r_4> > const & in, TMatrix< complex<r_4> > & out);
-  virtual void FFTBackward(TMatrix< complex<r_4> > const & in, TMatrix< complex<r_4> > & out);
-  virtual void FFTForward(TMatrix< r_4 > const & in, TMatrix< complex<r_4> > & out);
-  virtual void FFTBackward(TMatrix< complex<r_4> > const & in, TMatrix< r_4 > & out);
-
-
+  // Transforme N-dim sur des float
+  virtual void FFTForward(TArray< complex<r_4> > const & in, TArray< complex<r_4> > & out);
+  virtual void FFTBackward(TArray< complex<r_4> > const & in, TArray< complex<r_4> > & out);
+  virtual void FFTForward(TArray< r_4 > const & in, TArray< complex<r_4> > & out);
+  virtual void FFTBackward(TArray< complex<r_4> > const & in, TArray< r_4 > & out);
 
  protected:
@@ -66 +64 @@
 } // Fin du namespace
 
+template <class T> 
+class FFTArrayChecker {
+public:
+                FFTArrayChecker(bool checkpack=true, bool onedonly=false);
+  virtual       ~FFTArrayChecker();
+  virtual int   CheckResize(TArray< complex<T> > const & in, TArray< complex<T> > & out);
+  virtual int   CheckResize(TArray< T > const & in, TArray< complex<T> > & out);
+  virtual int   CheckResize(TArray< complex<T> > const & in, TArray< T > & out);
+  virtual void  ReShapetoReal(TArray< complex<T> > const & in, TArray< T >  & out);
+  virtual void  ReShapetoCompl(TArray< T > const & in, TArray< complex<T> > & out); 
+
+protected:
+  bool _checkpack;
+  bool _onedonly;
+};
+
 #endif