| 1 | #include "fftservintf.h" | 
|---|
| 2 |  | 
|---|
| 3 |  | 
|---|
| 4 | /*! | 
|---|
| 5 | \class SOPHYA::FFTServerInterface | 
|---|
| 6 | \ingroup NTools | 
|---|
| 7 | Defines the interface for FFT (Fast Fourier Transform) operations. | 
|---|
| 8 | */ | 
|---|
| 9 |  | 
|---|
| 10 | /* --Methode-- */ | 
|---|
| 11 | FFTServerInterface::FFTServerInterface(string info) | 
|---|
| 12 | { | 
|---|
| 13 | _info = info; | 
|---|
| 14 | _fgnorm = true; | 
|---|
| 15 | } | 
|---|
| 16 |  | 
|---|
| 17 | /* --Methode-- */ | 
|---|
| 18 | FFTServerInterface::~FFTServerInterface() | 
|---|
| 19 | { | 
|---|
| 20 | } | 
|---|
| 21 |  | 
|---|
| 22 | // ----------------- Transforme pour les double ------------------- | 
|---|
| 23 |  | 
|---|
| 24 | /* --Methode-- */ | 
|---|
| 25 | void FFTServerInterface::FFTForward(TArray< complex<r_8> > const &, TArray< complex<r_8> > &) | 
|---|
| 26 | { | 
|---|
| 27 | throw NotAvailableOperation("FFTServer::FFTForward(TArray...) Unsupported operation !"); | 
|---|
| 28 | } | 
|---|
| 29 |  | 
|---|
| 30 | /* --Methode-- */ | 
|---|
| 31 | void FFTServerInterface::FFTBackward(TArray< complex<r_8> > const &, TArray< complex<r_8> > &) | 
|---|
| 32 | { | 
|---|
| 33 | throw NotAvailableOperation("FFTServer::FFTBackward(TArray...) Unsupported operation !"); | 
|---|
| 34 | } | 
|---|
| 35 |  | 
|---|
| 36 | /* --Methode-- */ | 
|---|
| 37 | void FFTServerInterface::FFTForward(TArray< r_8 > const &, TArray< complex<r_8> > &) | 
|---|
| 38 | { | 
|---|
| 39 | throw NotAvailableOperation("FFTServer::FFTForward(TArray...) Unsupported operation !"); | 
|---|
| 40 | } | 
|---|
| 41 |  | 
|---|
| 42 | /* --Methode-- */ | 
|---|
| 43 | void FFTServerInterface::FFTBackward(TArray< complex<r_8> > const &, TArray< r_8 > &) | 
|---|
| 44 | { | 
|---|
| 45 | throw NotAvailableOperation("FFTServer::FFTBackward(TArray...) Unsupported operation !"); | 
|---|
| 46 | } | 
|---|
| 47 |  | 
|---|
| 48 |  | 
|---|
| 49 | // ----------------- Transforme pour les float ------------------- | 
|---|
| 50 |  | 
|---|
| 51 | /* --Methode-- */ | 
|---|
| 52 | void FFTServerInterface::FFTForward(TArray< complex<r_4> > const &, TArray< complex<r_4> > &) | 
|---|
| 53 | { | 
|---|
| 54 | throw NotAvailableOperation("FFTServer::FFTForward(TArray r_4 ... ) Unsupported operation !"); | 
|---|
| 55 | } | 
|---|
| 56 |  | 
|---|
| 57 | /* --Methode-- */ | 
|---|
| 58 | void FFTServerInterface::FFTBackward(TArray< complex<r_4> > const &, TArray< complex<r_4> > &) | 
|---|
| 59 | { | 
|---|
| 60 | throw NotAvailableOperation("FFTServer::FFTBackward(TArray r_4 ... ) Unsupported operation !"); | 
|---|
| 61 | } | 
|---|
| 62 |  | 
|---|
| 63 | /* --Methode-- */ | 
|---|
| 64 | void FFTServerInterface::FFTForward(TArray< r_4 > const &, TArray< complex<r_4> > &) | 
|---|
| 65 | { | 
|---|
| 66 | throw NotAvailableOperation("FFTServer::FFTForward(TArray r_4 ... ) Unsupported operation !"); | 
|---|
| 67 | } | 
|---|
| 68 |  | 
|---|
| 69 | /* --Methode-- */ | 
|---|
| 70 | void FFTServerInterface::FFTBackward(TArray< complex<r_4> > const &, TArray< r_4 > &) | 
|---|
| 71 | { | 
|---|
| 72 | throw NotAvailableOperation("FFTServer::FFTBackward(TArray r_4 ... ) Unsupported operation !"); | 
|---|
| 73 | } | 
|---|
| 74 |  | 
|---|
| 75 |  | 
|---|
| 76 |  | 
|---|
| 77 | /* --Methode-- */ | 
|---|
| 78 | template <class T> | 
|---|
| 79 | FFTArrayChecker<T>::FFTArrayChecker(string msg, bool checkpack, bool onedonly) | 
|---|
| 80 | { | 
|---|
| 81 | _msg = msg + " FFTArrayChecker::"; | 
|---|
| 82 | _checkpack = checkpack; | 
|---|
| 83 | _onedonly = onedonly; | 
|---|
| 84 | } | 
|---|
| 85 |  | 
|---|
| 86 | /* --Methode-- */ | 
|---|
| 87 | template <class T> | 
|---|
| 88 | FFTArrayChecker<T>::~FFTArrayChecker() | 
|---|
| 89 | { | 
|---|
| 90 | } | 
|---|
| 91 |  | 
|---|
| 92 | template <class T> | 
|---|
| 93 | T FFTArrayChecker<T>::ZeroThreshold() | 
|---|
| 94 | { | 
|---|
| 95 | return(0); | 
|---|
| 96 | } | 
|---|
| 97 |  | 
|---|
| 98 | r_8 FFTArrayChecker< r_8 >::ZeroThreshold() | 
|---|
| 99 | { | 
|---|
| 100 | return(1.e-18); | 
|---|
| 101 | } | 
|---|
| 102 |  | 
|---|
| 103 | r_4 FFTArrayChecker< r_4 >::ZeroThreshold() | 
|---|
| 104 | { | 
|---|
| 105 | return(1.e-9); | 
|---|
| 106 | } | 
|---|
| 107 |  | 
|---|
| 108 |  | 
|---|
| 109 |  | 
|---|
| 110 | /* --Methode-- */ | 
|---|
| 111 | template <class T> | 
|---|
| 112 | int FFTArrayChecker<T>::CheckResize(TArray< complex<T> > const & in, TArray< complex<T> > & out) | 
|---|
| 113 | { | 
|---|
| 114 | int k; | 
|---|
| 115 | string msg; | 
|---|
| 116 | if (in.Size() < 1) { | 
|---|
| 117 | msg = _msg + "CheckResize(complex in, complex out) - Unallocated input array !"; | 
|---|
| 118 | throw(SzMismatchError(msg)); | 
|---|
| 119 | } | 
|---|
| 120 | if (_checkpack) | 
|---|
| 121 | if ( !in.IsPacked() ) { | 
|---|
| 122 | msg = _msg + "CheckResize(complex in, complex out) - Not packed input array !"; | 
|---|
| 123 | throw(SzMismatchError(msg)); | 
|---|
| 124 | } | 
|---|
| 125 | int ndg1 = 0; | 
|---|
| 126 | for(k=0; k<in.NbDimensions(); k++) | 
|---|
| 127 | if (in.Size(k) > 1)  ndg1++; | 
|---|
| 128 | if (_onedonly) | 
|---|
| 129 | if (ndg1 > 1) { | 
|---|
| 130 | msg = _msg + "CheckResize(complex in, complex out) - Only 1-D array accepted !"; | 
|---|
| 131 | throw(SzMismatchError(msg)); | 
|---|
| 132 | } | 
|---|
| 133 | out.ReSize(in); | 
|---|
| 134 | //  sa_size_t sz[BASEARRAY_MAXNDIMS]; | 
|---|
| 135 | //  for(k=0; k<in.NbDimensions(); k++) | 
|---|
| 136 | //    sz[k] = in.Size(k); | 
|---|
| 137 | //  out.ReSize(in.NbDimensions(), sz); | 
|---|
| 138 |  | 
|---|
| 139 | return(ndg1); | 
|---|
| 140 | } | 
|---|
| 141 |  | 
|---|
| 142 | /* --Methode-- */ | 
|---|
| 143 | template <class T> | 
|---|
| 144 | int FFTArrayChecker<T>::CheckResize(TArray< T > const & in, TArray< complex<T> > & out) | 
|---|
| 145 | { | 
|---|
| 146 | int k; | 
|---|
| 147 | string msg; | 
|---|
| 148 | if (in.Size() < 1) { | 
|---|
| 149 | msg = _msg + "CheckResize(real in, complex out) - Unallocated input array !"; | 
|---|
| 150 | throw(SzMismatchError(msg)); | 
|---|
| 151 | } | 
|---|
| 152 | if (_checkpack) | 
|---|
| 153 | if ( !in.IsPacked() ) { | 
|---|
| 154 | msg = _msg + "CheckResize(real in, complex out) - Not packed input array !"; | 
|---|
| 155 | throw(SzMismatchError(msg)); | 
|---|
| 156 | } | 
|---|
| 157 | int ndg1 = 0; | 
|---|
| 158 | for(k=0; k<in.NbDimensions(); k++) | 
|---|
| 159 | if (in.Size(k) > 1)  ndg1++; | 
|---|
| 160 | if (_onedonly) | 
|---|
| 161 | if (ndg1 > 1) { | 
|---|
| 162 | msg = _msg + "CheckResize(real in, complex out) - Only 1-D array accepted !"; | 
|---|
| 163 | throw(SzMismatchError(msg)); | 
|---|
| 164 | } | 
|---|
| 165 | sa_size_t sz[BASEARRAY_MAXNDIMS]; | 
|---|
| 166 | for(k=0; k<in.NbDimensions(); k++) | 
|---|
| 167 | sz[k] = in.Size(k)/2+1; | 
|---|
| 168 | //    sz[k] = (in.Size(k)%2 != 0) ? in.Size(k)/2+1 : in.Size(k)/2; | 
|---|
| 169 |  | 
|---|
| 170 | out.ReSize(in.NbDimensions(), sz); | 
|---|
| 171 |  | 
|---|
| 172 | return(ndg1); | 
|---|
| 173 | } | 
|---|
| 174 |  | 
|---|
| 175 | /* --Methode-- */ | 
|---|
| 176 | template <class T> | 
|---|
| 177 | int FFTArrayChecker<T>::CheckResize(TArray< complex<T> > const & in, TArray< T > & out) | 
|---|
| 178 | { | 
|---|
| 179 | int k; | 
|---|
| 180 | string msg; | 
|---|
| 181 | if (in.Size() < 1) { | 
|---|
| 182 | msg = _msg + "CheckResize(complex in, real out) - Unallocated input array !"; | 
|---|
| 183 | throw(SzMismatchError(msg)); | 
|---|
| 184 | } | 
|---|
| 185 | if (_checkpack) | 
|---|
| 186 | if ( !in.IsPacked() ) { | 
|---|
| 187 | msg = _msg + "CheckResize(complex in, real out) - Not packed input array !"; | 
|---|
| 188 | throw(SzMismatchError(msg)); | 
|---|
| 189 | } | 
|---|
| 190 | int ndg1 = 0; | 
|---|
| 191 | for(k=0; k<in.NbDimensions(); k++) | 
|---|
| 192 | if (in.Size(k) > 1)  ndg1++; | 
|---|
| 193 | if (_onedonly) | 
|---|
| 194 | if (ndg1 > 1) { | 
|---|
| 195 | msg = _msg + "CheckResize(complex in, real out) - Only 1-D array accepted !"; | 
|---|
| 196 | throw(SzMismatchError(msg)); | 
|---|
| 197 | } | 
|---|
| 198 | sa_size_t sz[BASEARRAY_MAXNDIMS]; | 
|---|
| 199 | if (ndg1 > 1) { | 
|---|
| 200 | for(k=0; k<in.NbDimensions(); k++) | 
|---|
| 201 | sz[k] = in.Size(k)*2-1; | 
|---|
| 202 | } | 
|---|
| 203 | else { | 
|---|
| 204 | for(k=0; k<BASEARRAY_MAXNDIMS; k++)  sz[k] = 1; | 
|---|
| 205 | T thr = ZeroThreshold(); | 
|---|
| 206 | sa_size_t n = in.Size(in.MaxSizeKA()); | 
|---|
| 207 | sa_size_t ncs = ( (in[n-1].imag() < -thr) || (in[n-1].imag() > thr) ) ? | 
|---|
| 208 | ncs = 2*n-1 : ncs = 2*n-2; | 
|---|
| 209 | sz[in.MaxSizeKA()] = ncs; | 
|---|
| 210 | } | 
|---|
| 211 |  | 
|---|
| 212 | out.ReSize(in.NbDimensions(), sz); | 
|---|
| 213 |  | 
|---|
| 214 | return(ndg1); | 
|---|
| 215 |  | 
|---|
| 216 | } | 
|---|
| 217 |  | 
|---|
| 218 | /* --Methode-- | 
|---|
| 219 | template <class T> | 
|---|
| 220 | void FFTArrayChecker<T>::ReShapetoReal( TArray< complex<T> > const & in, TArray< T >  & out) | 
|---|
| 221 | { | 
|---|
| 222 | sa_size_t n = in.Size(); | 
|---|
| 223 | //  int ncs = (fabs(in(n-1).imag()) > 1.e-12) ? ncs = 2*n-1 : ncs = n*2-2; | 
|---|
| 224 | sa_size_t ncs = 2*n-1; | 
|---|
| 225 | sa_size_t k; | 
|---|
| 226 | out[0] = in[0].real(); | 
|---|
| 227 | for(k=1;k<n-1;k++) { | 
|---|
| 228 | out[2*k-1] = in[k].real(); | 
|---|
| 229 | out[2*k] = in[k].imag(); | 
|---|
| 230 | } | 
|---|
| 231 | //  if (ncs == n*2-2)  out[ncs-1] = in[n-1].real(); | 
|---|
| 232 | //  else { out[ncs-2] = in[n-1].real();  out[ncs-1] = in[n-1].imag(); } | 
|---|
| 233 | out[ncs-2] = in[n-1].real();  out[ncs-1] = in[n-1].imag(); | 
|---|
| 234 | } | 
|---|
| 235 | */ | 
|---|
| 236 |  | 
|---|
| 237 | /* --Methode-- | 
|---|
| 238 | template <class T> | 
|---|
| 239 | void FFTArrayChecker<T>::ReShapetoCompl(TArray< T > const & in, TArray< complex<T> > & out) | 
|---|
| 240 | { | 
|---|
| 241 | sa_size_t n =  in.Size(); | 
|---|
| 242 | sa_size_t ncs = n/2+1; | 
|---|
| 243 | sa_size_t nc = (n%2 != 0) ? n/2+1 : n/2; | 
|---|
| 244 | out[0] = complex<T> (in[0],0.); | 
|---|
| 245 | for(int k=1;k<nc;k++) | 
|---|
| 246 | out[k] =  complex<T> (in[2*k-1], in[2*k]); | 
|---|
| 247 | if (n%2 == 0) out[ncs-1] = complex<r_8>(in[n-1], 0.); | 
|---|
| 248 | } | 
|---|
| 249 | */ | 
|---|
| 250 |  | 
|---|
| 251 | #ifdef __CXX_PRAGMA_TEMPLATES__ | 
|---|
| 252 | #pragma define_template FFTArrayChecker<r_4> | 
|---|
| 253 | #pragma define_template FFTArrayChecker<r_8> | 
|---|
| 254 | #endif | 
|---|
| 255 |  | 
|---|
| 256 | #if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES) | 
|---|
| 257 | template class FFTArrayChecker<r_4>; | 
|---|
| 258 | template class FFTArrayChecker<r_8>; | 
|---|
| 259 | #endif | 
|---|