C:/temp/src/j2k/Deprecated/Complex.old/Complex.old Source File

00001 #ifndef __J2K__JComplex_CPP__
00002 #define __J2K__JComplex_CPP__
00003 
00004 #include <j2k/Fred/Math/Complex.hpp>
00005 
00006 /////////////////////////////////////////////////////////////////////////////
00007 /// Constructors and Destructors                                          ///
00008 /////////////////////////////////////////////////////////////////////////////
00009 inline JComplex::JComplex( const CplxDouble r, const CplxDouble i )
00010   : Real( r ), Imag( i ) { }
00011 
00012 inline JComplex::JComplex( const CplxDouble r )
00013   : Real( r ), Imag( 0.0 ) { }
00014 
00015 inline JComplex::JComplex( const JComplex& c )
00016   : Real( c.Real ), Imag( c.Imag ) { }
00017 
00018 inline JComplex::~JComplex() { }
00019 
00020 /////////////////////////////////////////////////////////////////////////////
00021 /// Compare                                                               ///
00022 /////////////////////////////////////////////////////////////////////////////
00023 static BOOL JComplex::Compare( const JComplex& c ) const {
00024   if (this == &c) return 1;
00025   return ( (Real == c.Real)  &&  (Imag == c.Imag) );
00026 }
00027 
00028 inline BOOL JComplex::operator==( const JComplex& c ) const {
00029   return Compare( c );
00030 }
00031 
00032 inline BOOL JComplex::operator!=( const JComplex& c ) const {
00033   return !Compare( c );
00034 }
00035 
00036 static BOOL JComplex::Compare( const CplxDouble r ) const {
00037   return ( (Real == r)  &&  (Imag == 0) );
00038 }
00039 
00040 inline BOOL JComplex::operator==( const CplxDouble r ) const {
00041   return Compare( r );
00042 }
00043 
00044 inline BOOL JComplex::operator!=( const CplxDouble r ) const {
00045   return !Compare( r );
00046 }
00047 
00048 static BOOL JComplex::Compare( const JComplex& l, const JComplex& r ) {
00049   return ( (l.Real == r.Real) && (l.Imag == r.Imag) );                   
00050 }                                                                        
00051                                                                            
00052 inline friend BOOL JComplex::operator==( const JComplex& l, const JComplex& r ) {  
00053   return Compare( l, r );                                                
00054 }                                                                        
00055                                                                            
00056 inline friend BOOL JComplex::operator!=( const Number& l, const Number& r ) {      
00057   return !Compare( l, r );                                               
00058 }                                                                        
00059 
00060 /////////////////////////////////////////////////////////////////////////////
00061 /// Utilities                                                             ///
00062 /////////////////////////////////////////////////////////////////////////////
00063 inline JComplex JComplex::getConjugate() const {
00064   return JComplex( Real, 0.0 - Imag );
00065 }
00066 
00067 inline CplxDouble JComplex::Real() const {
00068   return Real;
00069 }
00070 
00071 inline CplxDouble JComplex::Imag() const {
00072   return Imag;
00073 }
00074 
00075 inline CplxDouble JComplex::Real( CplxDouble x ) const {
00076   return (Real = x);
00077 }
00078 
00079 inline CplxDouble JComplex::Imag( CplxDouble y ) const {
00080   return (Imag = y);
00081 }
00082 
00083 inline CplxDouble JComplex::Modulus() const {
00084   return ((Real * Real) + (Imag * Imag));
00085 }
00086 
00087 inline void JComplex::Conjugate() {
00088   // Negate the imaginary component
00089   Imag = 0.0 - Imag;
00090 }
00091 
00092 /////////////////////////////////////////////////////////////////////////////
00093 /// Assignment, Addition, Substraction operations                         ///
00094 /////////////////////////////////////////////////////////////////////////////
00095 static JComplex& JComplex::operator=( const CplxNumber& x ) {                            
00096   Real = x;                                                            
00097   Imag = 0;                                                            
00098   return (*this);                                                      
00099 }                                                                       
00100                                                                            
00101 static JComplex& JComplex::operator=( const JComplex& c ) {
00102   Real = c.Real;                                                       
00103   Imag = c.Imag;                                                       
00104   return (*this);                                                      
00105 }                                                                       
00106 
00107 inline JComplex& JComplex::operator+=( const JComplex& c ) {
00108   Real = Real + c.Real;
00109   Imag = Imag + c.Imag;
00110   return (*this);
00111 }
00112 
00113 inline JComplex& JComplex::operator-=( const JComplex& c ) {
00114   Real = Real - c.Real;
00115   Imag = Imag - c.Imag;
00116   return (*this);
00117 }
00118 
00119 inline JComplex& operator+=( JComplex& x, const JComplex& y ) {
00120   x.Real = x.Real + y.Real;
00121   x.Imag = x.Imag + y.Imag;
00122   return (x);
00123 }
00124 
00125 inline JComplex& operator-=( JComplex& x, const JComplex& y ) {
00126   x.Real = x.Real - y.Real;
00127   x.Imag = x.Imag - y.Imag;
00128   return (x);
00129 }                                                                       
00130 
00131 inline JComplex& JComplex::operator+=( const CplxDouble r ) {
00132   Real = Real + r;
00133   return (*this);
00134 }
00135 
00136 inline JComplex& JComplex::operator-=( const CplxDouble r ) {
00137   Real = Real - r;
00138   return (*this);
00139 }
00140 
00141 inline JComplex operator+(const JComplex& left, const JComplex& right) {
00142   return JComplex( (left.Real + right.Real), (left.Imag + right.Imag) );
00143 }
00144 
00145 inline JComplex operator-(const JComplex& left, const JComplex& right) {
00146   return JComplex( (left.Real - right.Real), (left.Imag - right.Imag) );
00147 }
00148 
00149 /////////////////////////////////////////////////////////////////////////////
00150 /// Scaling operations                                                    ///
00151 /////////////////////////////////////////////////////////////////////////////
00152 inline void JComplex::operator*=(const CplxDouble r) {
00153   Real = Real * r;
00154   Imag = Imag * r;
00155 }
00156 
00157 inline void JComplex::operator/=(const CplxDouble r) {
00158   Real = Real / r;
00159   Imag = Imag / r;
00160 }
00161 
00162 /////////////////////////////////////////////////////////////////////////////
00163 /// Multiplication and Division operations                                ///
00164 /////////////////////////////////////////////////////////////////////////////
00165 JComplex& operator*=( JComplex& x, const JComplex& y ) {
00166    // Reminder: (a+bj)(c+dj) = (ac + bcj + adj + bdjj)
00167    // Result:   (ac - bd) + ( bc + ad )j
00168    // jj = -1, by definition.
00169 
00170    CplxNumber Temp_Re = (x.Real * y.Real) - (x.Imag * y.Imag);
00171    CplxNumber Temp_Im = (x.Imag * y.Real) + (x.Real * y.Imag);
00172 
00173    x.Real = Temp_Re;
00174    x.Imag = Temp_Im;
00175 
00176    return (x);
00177 }
00178 
00179 JComplex& operator/=( JComplex& x, const JComplex& y ) {
00180    // Reminder: (a+bj)/(c+dj) = (a+bj)(c-dj) / [ (c+dj)(c-dj) ]
00181    // = (ac + bcj - adj - bdjj) / (cc - ddjj)
00182    // = [ (ac + bd) + ( bc - ad )j ] / [ cc + dd ]
00183    // jj = -1, by definition.
00184 
00185    // Doesn't handle Divide by Zero, but STL does it.
00186 
00187    CplxNumber Temp_delta = (y.Real * y.Real) + (y.Imag * y.Imag);
00188    CplxNumber Temp_Re    = (x.Real * y.Real) + (x.Imag * y.Imag);
00189    CplxNumber Temp_Im    = (x.Imag * y.Real) - (x.Real * y.Imag);
00190 
00191    x.Real = Temp_Re / Temp_delta;
00192    x.Imag = Temp_Im / Temp_delta;
00193 
00194    return (x);
00195 }
00196 
00197 static void JComplex::operator*=(const JComplex& c) {
00198   // Don't affect either until both are done, do new real to temp
00199   const CplxDouble tmpReal = (Real * c.Real) - (Imag * c.Imag);
00200   Imag = (Real * c.Imag) + (Imag * c.Real);
00201   Real = tmpReal;
00202 }
00203 
00204 JComplex operator*(const JComplex& left, const JComplex& right) {
00205   return JComplex (
00206    (left.Real * right.Real) - (left.Imag * right.Imag),
00207    (left.Real * right.Imag) + (left.Imag * right.Real)
00208   );
00209 }
00210 
00211 void JComplex::operator/=(const JComplex& c) {
00212  const CplxDouble Tmp1 = c.Real / c.Imag;
00213  const CplxDouble Tmp2 = c.Imag + (Tmp1 * c.Real);
00214 
00215  // Don't affect either until both are done, do new real to temp
00216  const CplxDouble TmpReal = (Real * Tmp1 + Imag) / Tmp2;
00217  Imag = (Imag * Tmp1 - Real) / Tmp2;
00218  Real = TmpReal;
00219 }
00220 
00221 JComplex operator/(const JComplex& left, const JComplex& right) {
00222  const CplxDouble Tmp1 = right.Real / right.Imag;
00223  const CplxDouble Tmp2 = right.Imag + (Tmp1 * right.Real);
00224 
00225  return JComplex(
00226   ( (left.Real * Tmp1 + left.Imag) / Tmp2 ),
00227   ( (left.Imag * Tmp1 - left.Real) / Tmp2 )
00228  );
00229 }
00230 
00231 /////////////////////////////////////////////////////////////////////////////
00232 /// Utilities                                                             ///
00233 /////////////////////////////////////////////////////////////////////////////
00234    CplxNumber Real(const CplxNumber& x) {                                  
00235       return (Real = x);                                                   
00236    }                                                                       
00237                                                                            
00238    CplxNumber Imag(const CplxNumber& y) {                                  
00239       return (Imag = y);                                                   
00240    }                                                                       
00241                                                                            
00242    CplxNumber Real() const {                                               
00243       return (Real);                                                       
00244    }                                                                       
00245                                                                            
00246    CplxNumber Imag() const {                                               
00247       return (Imag);                                                       
00248    }                                                                       
00249                                                                            
00250    // Arithmetic Operations
00251    JComplex& operator=( const CplxNumber& x ) {                            
00252       Real = x;                                                            
00253       Imag = 0;                                                            
00254       return (*this);                                                      
00255    }                                                                       
00256                                                                            
00257    JComplex& operator=( const JComplex& c ) {
00258       Real = c.Real;                                                       
00259       Imag = c.Imag;                                                       
00260       return (*this);                                                      
00261    }                                                                       
00262                                                                            
00263    JComplex& operator+=(const CplxNumber& x ) {                            
00264       Real = Real + x;                                                     
00265       return (*this);                                                      
00266    }                                                                       
00267                                                                            
00268    JComplex& operator+=( JComplex& x, const JComplex& y ) {                
00269      x.Real = x.Real + y.Real;
00270      x.Imag = y.Imag + y.Imag;
00271      return (x);                                                           
00272    }                                                                       
00273                                                                            
00274    JComplex& operator-=(const CplxNumber& x) {                             
00275       Real = Real - x;                                                     
00276       return (*this);                                                      
00277    }                                                                       
00278                                                                            
00279    JComplex& operator-=( JComplex& x, const JComplex& y ) {                
00280      x.Real = x.Real - y.Real;
00281      x.Imag = y.Imag - y.Imag;
00282      return (x);                                                           
00283    }                                                                       
00284                                                                            
00285    JComplex& operator*=(const CplxNumber& x) {                             
00286       Real = Real * x;                                                     
00287       Imag = Imag * x;                                                     
00288       return (*this);                                                      
00289    }                                                                       
00290                                                                            
00291    JComplex& operator/=(const CplxNumber& x) {                             
00292       Real = Real / x;                                                     
00293       Imag = Imag / x;                                                     
00294       return (*this);                                                      
00295    }                                                                       
00296                                                                            
00297 protected:                                                                 
00298    CplxNumber Real;                                                        
00299    CplxNumber Imag;                                                        
00300                                                                            
00301 public:                                                                    
00302                                                                            
00303   static BOOL Compare( const JComplex& l, const JComplex& r ) {            
00304     return ( (l.Real == r.Real) && (l.Imag == r.Imag) );                   
00305   }                                                                        
00306                                                                            
00307   inline friend BOOL operator==( const JComplex& l, const JComplex& r ) {  
00308     return Compare( l, r );                                                
00309   }                                                                        
00310                                                                            
00311   inline friend BOOL operator!=( const Number& l, const Number& r ) {      
00312     return !Compare( l, r );                                               
00313   }                                                                        
00314                                                                            
00315                                                                            
00316 
00317    JComplex& operator*=( JComplex& x, const JComplex& y ) {
00318      // Reminder: (a+bj)(c+dj) = (ac + bcj + adj + bdjj)
00319      // Result:   (ac - bd) + ( bc + ad )j
00320      // jj = -1, by definition.
00321 
00322      CplxNumber Temp_Re = (x.Real * y.Real) - (x.Imag * y.Imag);
00323      CplxNumber Temp_Im = (x.Imag * y.Real) + (x.Real * y.Imag);
00324 
00325      x.Real = Temp_Re;
00326      x.Imag = Temp_Im;
00327 
00328      return (x);
00329    }
00330 
00331    JComplex& operator/=( JComplex& x, const JComplex& y ) {
00332      // Reminder: (a+bj)/(c+dj) = (a+bj)(c-dj) / [ (c+dj)(c-dj) ]
00333      // = (ac + bcj - adj - bdjj) / (cc - ddjj)
00334      // = [ (ac + bd) + ( bc - ad )j ] / [ cc + dd ]
00335      // jj = -1, by definition.
00336 
00337      // Doesn't handle Divide by Zero, but STL does it.
00338 
00339      CplxNumber Temp_delta = (y.Real * y.Real) + (y.Imag * y.Imag);
00340 
00341      CplxNumber Temp_Re = (x.Real * y.Real) + (x.Imag * y.Imag);
00342      CplxNumber Temp_Im = (x.Imag * y.Real) - (x.Real * y.Imag);
00343 
00344      x.Real = Temp_Re / Temp_delta;
00345      x.Imag = Temp_Im / Temp_delta;
00346 
00347      return (x);
00348    }
00349 
00350    static CplxNumber _Infv(CplxNumber) {                                   
00351      return ( _Inf._D );                                                   
00352    }                                                                       
00353                                                                            
00354    static bool _Isinf(CplxNumber x) {                                      
00355      double temp = (double)x;                                              
00356      return (_Dtest(&temp) == _INFCODE);                                   
00357    }                                                                       
00358                                                                            
00359    static bool _Isnan(CplxNumber x) {                                      
00360      double temp = (double)x;                                              
00361      return (_Dtest(&temp) == _NANCODE);                                   
00362    }                                                                       
00363                                                                            
00364    static CplxNumber _Nanv(CplxNumber) {                                   
00365      return (_Nan._D);                                                     
00366    }                                                                       
00367 
00368 /* STL content
00369 // It use Complex double as temp to verify for NaN, Inf
00370 // Might use Double for that.
00371 {
00372    if ( _Isnan( Temp_re ) || _Isnan( Temp_im ) ) {
00373      x.Real( _Nanv( Temp_re )), x.Imag( x.Real );
00374    } else if (
00375    (Temp_im < 0 ? -Temp_im : +Temp_im) < (Temp_re < 0 ? -Temp_re : +Temp_re)
00376    ) {
00377      CplxNumber _Wr = Temp_im / Temp_re;
00378      CplxNumber _Wd = Temp_re + _Wr * Temp_im;
00379 
00380      if (_Isnan(_Wd) || _Wd == 0) {
00381        x.Real(_Nanv(Temp_re)), x.Imag(x.Real());
00382      } else {
00383        CplxNumber temp = (x.Real() + x.Imag() * _Wr) / _Wd;
00384        x.Imag((x.Imag() - x.Real() * _Wr) / _Wd);
00385        x.Real(temp);
00386      }
00387 
00388    } else if (Temp_im == 0) {
00389       x.Real(_Nanv(Temp_re)), x.Imag(x.Real());
00390    } else {
00391       CplxNumber _Wr = Temp_re / Temp_im;
00392       CplxNumber _Wd = Temp_im + _Wr * Temp_re;
00393 
00394       if (_Isnan(_Wd) || _Wd == 0) {
00395          x.Real(_Nanv(Temp_re)), x.Imag(x.Real());
00396       } else {
00397          CplxNumber temp = (x.Real() * _Wr + x.Imag()) / _Wd;
00398          x.Imag((x.Imag() * _Wr - x.Real()) / _Wd);
00399          x.Real(temp);
00400       }
00401    }
00402    return (x);
00403 }
00404 */
00405 
00406 /*
00407 
00408 // TEMPLATE FUNCTION operator>>
00409 template<class Exponent, class _Tr, class _U> inline
00410    basic_istream<Exponent, _Tr>&  operator>>(
00411       basic_istream<Exponent, _Tr>& Im, JComplex& x)
00412 {
00413    Exponent _Ch;
00414    long double Real, Imag;
00415    if (Im >> _Ch && _Ch != '(')
00416       Im.putback(_Ch), Im >> Real, Imag = 0;
00417    else if (Im >> Real >> _Ch && _Ch != ',')
00418       if (_Ch == ')')
00419          Imag = 0;
00420       else
00421          Im.putback(_Ch), Im.setstate(ios_base::failbit);
00422    else if (Im >> Imag >> _Ch && _Ch != ')')
00423          Im.putback(_Ch), Im.setstate(ios_base::failbit);
00424    if (!Im.fail())
00425       x = JComplex((_U)Real, (_U)Imag);
00426    return (Im);
00427 }
00428 
00429 // TEMPLATE FUNCTION operator<<
00430 template<class Exponent, class _Tr, class _U> inline
00431    basic_ostream<Exponent, _Tr>&  operator<<(
00432       basic_ostream<Exponent, _Tr>& _O, const JComplex& x)
00433    {basic_ostringstream<Exponent, _Tr, allocator<Exponent> > _S;
00434    _S.flags(_O.flags());
00435    _S.imbue(_O.getloc());
00436    _S.precision(_O.precision());
00437    _S << '(' << Real(x) << ',' << Imag(x) << ')';
00438    return (_O << _S.str().c_str()); }
00439  #include <xcomplex>
00440 */
00441 
00442 #endif