00001
00002
00003
00004
00005
00006
00007 #include <iostream>
00008
00009 #ifdef STORM_USEF
00010 #include "rtmrartiwt_.h"
00011 #include "drtmrartiwt_.h"
00012 #include "irtmrartiwt_.h"
00013 #include "artmrartiwt_.h"
00014 #include "dartmrartiwt_.h"
00015 #include "iartmrartiwt_.h"
00016 #endif
00017
00018 #include "stormdef.h"
00019 #define INVSQRT2 STORMDEF_MATHCONST_INVSQRT2
00020
00021 void rtmrartiwt(
00022 const int& Msize,
00023 const int& mmax,
00024 int& mrastep,
00025 float* mracoef,
00026 float& data
00027 ) {
00028 #ifdef STORM_USEF
00029 rtmrartiwt_( &Msize, &mmax, &mrastep, mracoef, &data );
00030 #else
00031 mrastep++;
00032 int itmp = (1<<mmax) - 1;
00033 itmp = itmp - (itmp & mrastep);
00034 int imcycle=0;
00035 while ( ( itmp & (1<<imcycle) ) && ( imcycle < mmax ) ) { imcycle++; }
00036 int m = imcycle;
00037 if ( m < mmax ) {
00038 mracoef[2*m] = ( mracoef[2*(m+1)] - mracoef[2*(m+1)+1] ) * INVSQRT2;
00039 }
00040 m--;
00041 while ( m >= 0 ) {
00042 mracoef[2*m] = ( mracoef[2*(m+1)] + mracoef[2*(m+1)+1]) * INVSQRT2;
00043 m--;
00044 }
00045 data = mracoef[0];
00046 #endif
00047 }
00048
00049
00050 void rtmrartiwt(
00051 const int& Msize,
00052 const int& mmax,
00053 int& mrastep,
00054 double* mracoef,
00055 double& data
00056 ) {
00057 #ifdef STORM_USEF
00058 drtmrartiwt_( &Msize, &mmax, &mrastep, mracoef, &data );
00059 #else
00060 mrastep++;
00061 int itmp = (1<<mmax) - 1;
00062 itmp = itmp - ( itmp & mrastep);
00063 int imcycle = 0;
00064 while ( ( itmp & (1<<imcycle) ) && ( imcycle < mmax ) ) { imcycle++; }
00065 int m = imcycle;
00066 if ( m < mmax ) {
00067 mracoef[2*m] = ( mracoef[2*(m+1)] - mracoef[2*(m+1)+1]) * INVSQRT2;
00068 }
00069 m--;
00070 while ( m >= 0 ) {
00071 mracoef[2*m] = ( mracoef[2*(m+1)] + mracoef[2*(m+1)+1]) * INVSQRT2;
00072 m--;
00073 }
00074 data = mracoef[0];
00075 #endif
00076 }
00077
00078
00079 void rtmrartiwt(
00080 const int& Msize,
00081 const int& mmax,
00082 int& mrastep,
00083 int* mracoef,
00084 int& data
00085 ) {
00086 #ifdef STORM_USEF
00087 irtmrartiwt_( &Msize, &mmax, &mrastep, mracoef, &data );
00088 #else
00089 mrastep++;
00090 int itmp = (1<<mmax) - 1;
00091 itmp = itmp - (itmp & mrastep);
00092 int imcycle = 0;
00093 while ( ( itmp & (1<<imcycle) ) && ( imcycle < mmax ) ) { imcycle++; }
00094 int m = imcycle;
00095 if ( m < mmax ) {
00096 itmp = mracoef[2*(m+1)+1];
00097 if ( (itmp >= 0) || ((itmp & 0x00000001) == 0) ) {
00098 mracoef[2*m] = mracoef[2*(m+1)] - itmp/2;
00099 } else {
00100 mracoef[2*m] = mracoef[2*(m+1)] - (itmp/2-1);
00101 }
00102 }
00103 m--;
00104 while ( m >= 0 ) {
00105 itmp = mracoef[2*(m+1)+1] + 1;
00106 if ( (itmp >= 0) || ((itmp & 0x00000001) == 0) ) {
00107 mracoef[2*m] = mracoef[2*(m+1)] + itmp/2;
00108 } else {
00109 mracoef[2*m] = mracoef[2*(m+1)] + (itmp/2-1);
00110 }
00111 m--;
00112 }
00113 data = mracoef[0];
00114 #endif
00115 }
00116
00117
00118 void rtmrartiwt(
00119 const int& Asize,
00120 const int& Msize,
00121 const int& amax,
00122 const int& mmax,
00123 int& mrastep,
00124 float* mracoef,
00125 float* data
00126 ) {
00127 #ifdef STORM_USEF
00128 artmrartiwt_( &Asize,&Msize, &amax,&mmax, &mrastep, mracoef, data );
00129 #else
00130 mrastep++;
00131 int itmp = (1<<mmax) - 1;
00132 itmp = itmp - ( itmp & mrastep );
00133 int imcycle = 0;
00134 while ( ( itmp & (1 <<imcycle) ) && ( imcycle < mmax ) ) { imcycle++; }
00135 int m;
00136 for(int a=0; a <= amax; a++ ) {
00137 m=imcycle;
00138 if ( m < mmax ) {
00139 mracoef[ a*2*Msize + 2*m ]=( mracoef[ a*2*Msize + 2*(m+1) ]
00140 - mracoef[ a*2*Msize + 2*(m+1) + 1 ] ) * INVSQRT2;
00141 }
00142 m--;
00143 while ( m >= 0 ) {
00144 mracoef[ a*2*Msize + 2*m ] = ( mracoef[ a*2*Msize + 2*(m+1) ]
00145 +mracoef[ a*2*Msize + 2*(m+1) + 1 ] ) * INVSQRT2;
00146 m--;
00147 }
00148 data[a] = mracoef[a*2*Msize];
00149 }
00150 #endif
00151 }
00152
00153
00154 void rtmrartiwt(
00155 const int& Asize,
00156 const int& Msize,
00157 const int& amax,
00158 const int& mmax,
00159 int& mrastep,
00160 double* mracoef,
00161 double* data
00162 ) {
00163 #ifdef STORM_USEF
00164 dartmrartiwt_( &Asize,&Msize, &amax,&mmax, &mrastep, mracoef, data );
00165 #else
00166 mrastep++;
00167 int itmp = (1<<mmax) - 1;
00168 itmp = itmp - ( itmp & mrastep );
00169 int imcycle=0;
00170 while ( ( itmp & (1<<imcycle) ) && ( imcycle < mmax ) ) { imcycle++; }
00171 int m;
00172 for(int a = 0; a <= amax; a++ ) {
00173 m=imcycle;
00174 if ( m < mmax ) {
00175 mracoef[ a*2*Msize + 2*m ] = ( mracoef[ a*2*Msize + 2*(m+1) ]
00176 - mracoef[ a*2*Msize + 2*(m+1) + 1 ] ) * INVSQRT2;
00177 }
00178 m--;
00179 while ( m >= 0 ) {
00180 mracoef[ a*2*Msize + 2*m ] = ( mracoef[ a*2*Msize + 2*(m+1) ]
00181 +mracoef[ a*2*Msize + 2*(m+1) + 1 ] ) * INVSQRT2;
00182 m--;
00183 }
00184 data[a] = mracoef[a*2*Msize];
00185 }
00186 #endif
00187 }
00188
00189
00190 void rtmrartiwt(
00191 const int& Asize,
00192 const int& Msize,
00193 const int& amax,
00194 const int& mmax,
00195 int& mrastep,
00196 int* mracoef,
00197 int* data
00198 ) {
00199 #ifdef STORM_USEF
00200 iartmrartiwt_( &Asize,&Msize, &amax,&mmax, &mrastep, mracoef, data );
00201 #else
00202 mrastep++;
00203 int itmp = (1<<mmax) - 1;
00204 itmp = itmp - ( itmp & mrastep );
00205 int imcycle=0;
00206 while ( ( itmp & (1<<imcycle) ) && ( imcycle < mmax ) ) { imcycle++; }
00207 int m;
00208 for(int a = 0; a <= amax; a++ ) {
00209 m=imcycle;
00210 if ( m < mmax ) {
00211 itmp = mracoef[ a*2*Msize + 2*(m+1) + 1 ];
00212 if ( (itmp >= 0) || ((itmp & 0x00000001) == 0) ) {
00213 mracoef[a*2*Msize + 2*m]=mracoef[a*2*Msize + 2*(m+1)] - itmp/2;
00214 } else {
00215 mracoef[a*2*Msize + 2*m]=mracoef[a*2*Msize + 2*(m+1)] - (itmp/2-1);
00216 }
00217 }
00218 m--;
00219 while ( m >= 0 ) {
00220 itmp = mracoef[ a*2*Msize + 2*(m+1) + 1 ] + 1;
00221 if ( (itmp >= 0) || ((itmp & 0x00000001) == 0) ) {
00222 mracoef[a*2*Msize + 2*m]=mracoef[a*2*Msize + 2*(m+1)] + itmp/2;
00223 } else {
00224 mracoef[a*2*Msize + 2*m]=mracoef[a*2*Msize + 2*(m+1)] + (itmp/2-1);
00225 }
00226 m--;
00227 }
00228 data[a] = mracoef[a*2*Msize];
00229 }
00230 #endif
00231 }
00232
