C++ main module for emicrom Package  1.0
MATSGN_ComplexArray.h
Go to the documentation of this file.
1 #ifndef MATSGN_ComplexArray_H
2 #define MATSGN_ComplexArray_H
3 
4 #include "CORE_Object.h"
5 #include "CORE_Array.h"
6 #include "CORE_Exception.h"
7 #include "fftw3.h"
8 
16 class MATSGN_ComplexArray : public virtual CORE_Object {
18  // ATTRIBUTES
19 
20 
21 private:
22 
23  //values of the warray
24  fftw_complex *mValues;
25 
26  //size of the array
28 
29  //memory size of the array
31 
32  //inidate if the values has been allocated within this class
34 
35 
36  mutable tDouble mReal,mImag;
37 
38  // ASSOCIATIONS
39 
40 
41 public:
42  // METHODS
43 
44  // CONSTRUCTORS
45 
47  MATSGN_ComplexArray(void);
48 
50  MATSGN_ComplexArray(const tUIndex& n);
51 
52 
53 
54  // DESTRUCTORS
55 
56 
59  virtual ~MATSGN_ComplexArray(void);
60 
61 public:
62 
63  //operators
64  //=========
65 
70  inline fftw_complex& operator[](const tUIndex & i) {
71  ASSERT_IN(i<mSize);
72  return mValues[i];
73  };
78  inline const fftw_complex& operator[](const tUIndex & i) const {
79  ASSERT_IN(i<mSize);
80  return mValues[i];
81  };
82 
91  inline tDouble& operator()(const tUIndex & i,const tUSInt& j) {
92  ASSERT_IN(i<mSize);
93  ASSERT_IN(j<2);
94  return mValues[i][j];
95  };
103  inline const tDouble& operator()(const tUIndex & i,const tUSInt& j) const {
104  ASSERT_IN(i<mSize);
105  ASSERT_IN(j<2);
106  return mValues[i][j];
107  };
113  //cout << " complex array = operator \n";
114  copy(s);
115  return (*this);
116  };
121  //cout << " complex array = operator \n";
122  copy(s);
123  return (*this);
124  };
125 
131  copy(s);
132  return (*this);
133  };
139  copy(s);
140  return (*this);
141  };
146  inline MATSGN_ComplexArray& operator=(const tReal& s) {
147  initArray(s);
148  return (*this);
149  };
154  inline MATSGN_ComplexArray& operator=(const fftw_complex& s) {
155  initArray(s);
156  return (*this);
157  };
158 
164  if (mSize==0) return (*this);
165  tUIndex i;
166  fftw_complex* Ts=&(*this)[0];
167  for (i=0;i<mSize;i++) {
168  (*Ts)[0]*=s;
169  (*Ts)[1]*=s;
170  Ts++;
171 
172  }
173  return (*this);
174  };
181  if (mSize==0) return (*this);
182  if (Y.getSize()<mSize) throw CORE_Exception("math/signal",
183  "MATSGN_ComplexArray*=Y",
184  "Y has too small size "+CORE_Integer::toString(Y.getSize())+" with respect to This size :"+CORE_Integer::toString(mSize));
185  tUIndex i;
186  fftw_complex* Ts=&(*this)[0];
187  const fftw_complex *Ys=&Y[0];
188  for (i=0;i<mSize;i++) {
189  mReal=(*Ts)[0]*(*Ys)[0]-(*Ts)[1]*(*Ys)[1];
190  (*Ts)[1]=(*Ts)[0]*(*Ys)[1]+(*Ts)[1]*(*Ys)[0];
191  (*Ts)[0]=mReal;
192  Ts++;
193  Ys++;
194  }
195  return (*this);
196  };
197 
201  if (s==0) throw CORE_Exception("math/signam","MATSGN_ComplexArray/=","division by 0");
202  if (mSize==0) return (*this);
203  tUIndex i;
204  fftw_complex* Ts=&(*this)[0];
205  for (i=0;i<mSize;i++) {
206  (*Ts)[0]/=s;
207  (*Ts)[1]/=s;
208  Ts++;
209 
210  }
211  return (*this);
212  };
213 
219  if (mSize==0) return (*this);
220  tUIndex i;
221  fftw_complex* Ts=&(*this)[0];
222  for (i=0;i<mSize;i++) {
223  (*Ts)[0]+=s;
224  Ts++;
225 
226  }
227  return (*this);
228  };
234  if (mSize==0) return (*this);
235  tUIndex i;
236  fftw_complex* Ts=&(*this)[0];
237  for (i=0;i<mSize;i++) {
238  (*Ts)[0]-=s;
239  Ts++;
240 
241  }
242  return (*this);
243  };
244 
245  //Builder of array
246  //================
247 private:
248 
251  void desallocate();
252 
257  void allocate(const tUIndex& cap);
258 
259 public:
260  //copy method
261  //===========
262 
266  void copy(const MATSGN_ComplexArray& f);
267 
271  void copy(const CORE_ComplexArray& f);
272 
277  template<class T>
278  void copy(const complex<T>* f,const tUIndex& n) {
279  tUIndex k;
280 
281  if (n==0) {
282  desallocate();
283  return;
284  }
285 
286  //allocate the vector
287  if (n>mCapacity) {
288  mSize=n;
289  allocate(n);
290  }
291 
292  //copy the vector
293  if (n==0) return;
294 
295  const complex<T>* pf=&f[0];
296  fftw_complex* vs=mValues;
297  for (k=0;k<n;k++) {
298  (*vs)[0]=pf->real();
299  (*vs)[1]=pf->imag();
300  vs++;
301  pf++;
302 
303  }
304  }
305 
309  void copy(const CORE_RealArray& f);
310 
314  void copy(const CORE_DoubleArray& f);
315 
316 
317  // SET methods
318 
322  void initArray(const tReal& f);
323 
327  void initArray(const fftw_complex& f);
328 
329 
333  inline void setSize(const tUIndex& n) {
334  if (n>mCapacity) allocate(n);
335  mSize=n;
336  }
337 
340  void fitToSize();
341 
347  inline void setValue(const tUIndex& index,const tDouble& x,const tDouble& y) {
348  fftw_complex &c=mValues[index];
349  c[0]=x;
350  c[1]=y;
351  }
352 
359  inline tBoolean setValuesByReference(const tUIndex& n,fftw_complex* array) {
360  desallocate();
361  mValues=array;
362  mSize=n;
363  mCapacity=n;
364  mIsPrivateAllocated=false;
365  return true;
366  }
367 
368  // GET methods
369 
373  inline const fftw_complex* getValues() const {
374  return mValues;
375  }
379  inline fftw_complex* getValues() {
380  return mValues;
381  }
385  inline const tUIndex& getSize() const {
386  return mSize;
387  }
391  inline const tUIndex& getCapacity() const {
392  return mCapacity;
393  }
394 
395 
401  inline void swap(const tUIndex& i,const tUIndex& j) {
402  std::swap(mValues[i][0],mValues[j][0]);
403  std::swap(mValues[i][1],mValues[j][1]);
404  }
405 
410  tReal norm2() const;
411 
416  inline tReal norm() const {
417  return sqrt(norm2());
418  }
419 
420 
426  tReal distance2(const MATSGN_ComplexArray& x) const;
427 
428 
434  inline tReal distance(const MATSGN_ComplexArray& x) const {
435  return sqrt(distance2(x));
436  }
437 
438 
443  tReal fabs(tUIndex& i) const;
444 
450  tReal fabs(const MATSGN_ComplexArray& x,tUIndex& i) const;
451 
459  static void InPlaceTransposeWithMarker(fftw_complex* A, const tUIndex& nRows, const tUIndex& nCols);
460 
468  static void InPlaceTranspose(fftw_complex* A, const tUIndex& nRows, const tUIndex& nCols);
469 
472  virtual tString toString() const;
473 
474 
475 
476 
477 
478 };
479 
480 #endif
void copy(const complex< T > *f, const tUIndex &n)
copy the complex array
Definition: MATSGN_ComplexArray.h:278
#define tDouble
Definition: types.h:52
tReal distance(const MATSGN_ComplexArray &x) const
compute the distance of this to x
Definition: MATSGN_ComplexArray.h:434
MATSGN_ComplexArray & operator*=(const MATSGN_ComplexArray &Y)
multiply the array by Y
Definition: MATSGN_ComplexArray.h:180
virtual ~MATSGN_ComplexArray(void)
destroy an FFT complex array
Definition: MATSGN_ComplexArray.cpp:24
tReal distance2(const MATSGN_ComplexArray &x) const
compute the square of distance of this to x
Definition: MATSGN_ComplexArray.cpp:261
DEFINE_SPTR(MATSGN_ComplexArray)
tBoolean setValuesByReference(const tUIndex &n, fftw_complex *array)
set the the values by reference
Definition: MATSGN_ComplexArray.h:359
fftw_complex * getValues()
get values for reading or writing
Definition: MATSGN_ComplexArray.h:379
This class describes FFT complex array based on fft_complex structure.
Definition: MATSGN_ComplexArray.h:16
MATSGN_ComplexArray & operator=(const CORE_RealArray &s)
copy only the real part to s
Definition: MATSGN_ComplexArray.h:130
const fftw_complex & operator[](const tUIndex &i) const
get the value for reading only
Definition: MATSGN_ComplexArray.h:78
#define tUSInt
Definition: types.h:28
MATSGN_ComplexArray(void)
create a FFT complex array
Definition: MATSGN_ComplexArray.cpp:7
#define tBoolean
Definition: types.h:139
tString toString() const
return the string associated to the integer
Definition: CORE_Integer.h:106
const fftw_complex * getValues() const
get values for reading only
Definition: MATSGN_ComplexArray.h:373
static void InPlaceTranspose(fftw_complex *A, const tUIndex &nRows, const tUIndex &nCols)
transpose the vector A considered as a matrix of size nRows x nCols
Definition: MATSGN_ComplexArray.cpp:391
const tUIndex & getSize() const
get the size
Definition: MATSGN_ComplexArray.h:385
void swap(const tUIndex &i, const tUIndex &j)
swap the lines i & j
Definition: MATSGN_ComplexArray.h:401
MATSGN_ComplexArray & operator/=(const tReal &s)
divide operator
Definition: MATSGN_ComplexArray.h:200
void copy(const MATSGN_ComplexArray &f)
copy the complex array
Definition: MATSGN_ComplexArray.cpp:165
this class describes the exceptions raised for CORE package
Definition: CORE_Exception.h:15
this class describes an array
Definition: CORE_Array.h:19
MATSGN_ComplexArray & operator+=(const tReal &s)
add s to the array
Definition: MATSGN_ComplexArray.h:218
void desallocate()
desallocate the memory
Definition: MATSGN_ComplexArray.cpp:28
tBoolean mIsPrivateAllocated
Definition: MATSGN_ComplexArray.h:33
void initArray(const tReal &f)
init uniform value of f
Definition: MATSGN_ComplexArray.cpp:103
SP_OBJECT(MATSGN_ComplexArray)
MATSGN_ComplexArray & operator=(const tReal &s)
init uniform value
Definition: MATSGN_ComplexArray.h:146
void setSize(const tUIndex &n)
set size of the array
Definition: MATSGN_ComplexArray.h:333
tReal norm2() const
compute the square of norm
Definition: MATSGN_ComplexArray.cpp:231
tUIndex mSize
Definition: MATSGN_ComplexArray.h:27
#define tUIndex
Definition: types.h:126
static void InPlaceTransposeWithMarker(fftw_complex *A, const tUIndex &nRows, const tUIndex &nCols)
transpose the vector A considered as a matrix of size nRows x nCols
Definition: MATSGN_ComplexArray.cpp:367
abstract base class for most classes.
Definition: CORE_Object.h:53
#define tString
Definition: types.h:135
tUIndex mCapacity
Definition: MATSGN_ComplexArray.h:30
tReal norm() const
compute the norm
Definition: MATSGN_ComplexArray.h:416
MATSGN_ComplexArray & operator-=(const tReal &s)
sub s to the array
Definition: MATSGN_ComplexArray.h:233
void fitToSize()
set the capacityto size
Definition: MATSGN_ComplexArray.cpp:157
fftw_complex & operator[](const tUIndex &i)
get the value for reading & writing
Definition: MATSGN_ComplexArray.h:70
fftw_complex * mValues
Definition: MATSGN_ComplexArray.h:24
const tUIndex & getCapacity() const
get the capacity
Definition: MATSGN_ComplexArray.h:391
MATSGN_ComplexArray & operator=(const MATSGN_ComplexArray &s)
copy the values of the complex array
Definition: MATSGN_ComplexArray.h:112
MATSGN_ComplexArray & operator=(const fftw_complex &s)
init uniform value
Definition: MATSGN_ComplexArray.h:154
void allocate(const tUIndex &cap)
allocate the memory if the array
Definition: MATSGN_ComplexArray.cpp:43
MATSGN_ComplexArray & operator*=(const tReal &s)
multiply the array by s
Definition: MATSGN_ComplexArray.h:163
tReal fabs(tUIndex &i) const
compute the index of the array where the module is max
Definition: MATSGN_ComplexArray.cpp:244
tDouble mImag
Definition: MATSGN_ComplexArray.h:36
tDouble mReal
Definition: MATSGN_ComplexArray.h:36
const tDouble & operator()(const tUIndex &i, const tUSInt &j) const
get the value for reading only
Definition: MATSGN_ComplexArray.h:103
virtual tString toString() const
to string
Definition: MATSGN_ComplexArray.cpp:345
tDouble & operator()(const tUIndex &i, const tUSInt &j)
get the value for reading & writing
Definition: MATSGN_ComplexArray.h:91
void setValue(const tUIndex &index, const tDouble &x, const tDouble &y)
set the value at index
Definition: MATSGN_ComplexArray.h:347
#define tReal
Definition: types.h:118
MATSGN_ComplexArray & operator=(const CORE_ComplexArray &s)
copy the values of the tComplex array
Definition: MATSGN_ComplexArray.h:120
#define ASSERT_IN(a)
Definition: types.h:196
MATSGN_ComplexArray & operator=(const CORE_DoubleArray &s)
copy only the real part to s
Definition: MATSGN_ComplexArray.h:138