00001
00002
00003
00004
00005
00006
00007 #ifdef STORM_USEF
00008 #include "pd4fwt2ds_.h"
00009 #include "dpd4fwt2ds_.h"
00010 #include "pd4fwt2dsws_.h"
00011 #include "dpd4fwt2dsws_.h"
00012 #endif
00013
00014 #include "stormdef.h"
00015 #define H0 STORMDEF_MATHCONST_D4H0
00016 #define H1 STORMDEF_MATHCONST_D4H1
00017 #define H2 STORMDEF_MATHCONST_D4H2
00018 #define H3 STORMDEF_MATHCONST_D4H3
00019
00020 int pd4fwt2ds(
00021 const int& Sizef2d,
00022 const int& ns,
00023 const int& ks,
00024 const int& size,
00025 const float* sdata,
00026 float* wdata
00027 ) {
00028
00029 #ifdef STORM_USEF
00030
00031 int iflag;
00032 pd4fwt2ds_( &Sizef2d, &ns, &ks, &size, sdata, wdata, &iflag );
00033 return iflag;
00034
00035 #else
00036
00037 int ii,l;
00038 int l1,l2,l3;
00039 int size2;
00040 int size1 = size;
00041 float *u = new float[Sizef2d];
00042 float *v = new float[Sizef2d];
00043
00044 #ifdef STORM_FWT_CCHECKUSAGE
00045 if ( size1 < ns ) return 1;
00046 #endif
00047
00048 for( int i = 0; i < size1; i++ ) {
00049 ii = i * size;
00050 for( int j = 0; j < size1; j++ ) {
00051 wdata[ii+j] = sdata[ii+j];
00052 }
00053 }
00054
00055 if ( size1 == ns ) return 0;
00056
00057 while ( size1 > ns ) {
00058
00059 size2 = size1 / 2;
00060
00061 for( int j = 0; j < size1; j++ ) {
00062 for( int i = 0; i < size1; i++ ) {
00063 ii = i * size;
00064 u[i] = wdata[ii+j];
00065 }
00066 for( int i = 0; i < size2; i++ ) {
00067 l = 2 * i;
00068 l1 = (l+1)%size1;
00069 l2 = (l+2)%size1;
00070 l3 = (l+3)%size1;
00071 v[i] = H0*u[l] + H1*u[l1] + H2*u[l2] + H3*u[l3];
00072 v[size2+i] = H3*u[l] - H2*u[l1] + H1*u[l2] - H0*u[l3];
00073 }
00074 for( int i = 0; i < size1; i++ ) {
00075 ii = i * size;
00076 wdata[ii+j] = v[i];
00077 }
00078 }
00079
00080 for( int i = 0; i < size1; i++ ) {
00081 ii = i * size;
00082 for( int j = 0; j < size1; j++ ) {
00083 u[j] = wdata[ii+j];
00084 }
00085 for( int j = 0; j < size2; j++ ) {
00086 l = 2 * j;
00087 l1 = (l+1)%size1;
00088 l2 = (l+2)%size1;
00089 l3 = (l+3)%size1;
00090 v[j] = H0*u[l] + H1*u[l1] + H2*u[l2] + H3*u[l3];
00091 v[size2+j] = H3*u[l] - H2*u[l1] + H1*u[l2] - H0*u[l3];
00092 }
00093 for( int j = 0; j < size1; j++ ) {
00094 wdata[ii+j] = v[j];
00095 }
00096 }
00097
00098 size1 = size2;
00099
00100 }
00101
00102 delete [] u;
00103 delete [] v;
00104
00105 #ifdef STORM_FWT_CCHECKUSAGE
00106 if (size1 != ns) return 2;
00107 #endif
00108
00109 return 0;
00110
00111 #endif
00112 }
00113
00114
00115 int pd4fwt2ds(
00116 const int& Sizef2d,
00117 const int& ns,
00118 const int& ks,
00119 const int& size,
00120 const double* sdata,
00121 double* wdata
00122 ) {
00123
00124 #ifdef STORM_USEF
00125
00126 int iflag;
00127 dpd4fwt2ds_( &Sizef2d, &ns, &ks, &size, sdata, wdata, &iflag );
00128 return iflag;
00129
00130 #else
00131
00132 int ii,jj,l;
00133 int l1,l2,l3;
00134 int size2;
00135 int size1 = size;
00136 double *u = new double[Sizef2d];
00137 double *v = new double[Sizef2d];
00138
00139 #ifdef STORM_FWT_CCHECKUSAGE
00140 if ( size1 < ns ) return 1;
00141 #endif
00142
00143 for( int i = 0; i < size1; i++ ) {
00144 ii = i * size;
00145 for( int j = 0; j < size1; j++ ) {
00146 wdata[ii+j] = sdata[ii+j];
00147 }
00148 }
00149
00150 if ( size1 == ns ) return 0;
00151
00152 while ( size1 > ns ) {
00153
00154 size2 = size1 / 2;
00155
00156 for( int j = 0; j < size1; j++ ) {
00157 for( int i = 0; i < size1; i++ ) {
00158 ii = i * size;
00159 u[i] = wdata[ii+j];
00160 }
00161 for( int i = 0; i < size2; i++ ) {
00162 l = 2 * i;
00163 l1 = (l+1)%size1;
00164 l2 = (l+2)%size1;
00165 l3 = (l+3)%size1;
00166 v[i] = H0*u[l] + H1*u[l1] + H2*u[l2] + H3*u[l3];
00167 v[size2+i] = H3*u[l] - H2*u[l1] + H1*u[l2] - H0*u[l3];
00168 }
00169 for( int i = 0; i < size1; i++ ) {
00170 ii = i * size;
00171 wdata[ii+j] = v[i];
00172 }
00173 }
00174
00175 for( int i = 0; i < size1; i++ ) {
00176 ii = i * size;
00177 for( int j = 0; j < size1; j++ ) {
00178 u[j] = wdata[ii+j];
00179 }
00180 for( int j = 0; j < size2; j++ ) {
00181 l = 2 * j;
00182 l1 = (l+1)%size1;
00183 l2 = (l+2)%size1;
00184 l3 = (l+3)%size1;
00185 v[j] = H0*u[l] + H1*u[l1] + H2*u[l2] + H3*u[l3];
00186 v[size2+j] = H3*u[l] - H2*u[l1] + H1*u[l2] - H0*u[l3];
00187 }
00188 for( int j = 0; j < size1; j++ ) {
00189 jj = j * size;
00190 wdata[ii+j] = v[j];
00191 }
00192 }
00193
00194 size1 = size2;
00195
00196 }
00197
00198 delete [] u;
00199 delete [] v;
00200
00201 #ifdef STORM_FWT_CCHECKUSAGE
00202 if (size1 != ns) return 2;
00203 #endif
00204
00205 return 0;
00206
00207 #endif
00208 }
00209
00210
00211
00212 int pd4fwt2ds(
00213 const int& Sizef2d,
00214 const int& ns,
00215 const int& ks,
00216 const int& size,
00217 const float* sdata,
00218 float* wdata,
00219 float* ws
00220 ) {
00221
00222 #ifdef STORM_USEF
00223
00224 int iflag;
00225 pd4fwt2dsws_( &Sizef2d, &ns, &ks, &size, sdata, wdata, ws, &iflag );
00226 return iflag;
00227
00228 #else
00229
00230 int ii,l;
00231 int l1,l2,l3;
00232 int size2;
00233 int size1 = size;
00234 float *u = ws;
00235 float *v = &ws[Sizef2d];
00236
00237 #ifdef STORM_FWT_CCHECKUSAGE
00238 if ( size1 < ns ) return 1;
00239 #endif
00240
00241 for( int i = 0; i < size1; i++ ) {
00242 ii = i * size;
00243 for( int j = 0; j < size1; j++ ) {
00244 wdata[ii+j] = sdata[ii+j];
00245 }
00246 }
00247
00248 if ( size1 == ns ) return 0;
00249
00250 while ( size1 > ns ) {
00251
00252 size2 = size1 / 2;
00253
00254 for( int j = 0; j < size1; j++ ) {
00255 for( int i = 0; i < size1; i++ ) {
00256 ii = i * size;
00257 u[i] = wdata[ii+j];
00258 }
00259 for( int i = 0; i < size2; i++ ) {
00260 l = 2 * i;
00261 l1 = (l+1)%size1;
00262 l2 = (l+2)%size1;
00263 l3 = (l+3)%size1;
00264 v[i] = H0*u[l] + H1*u[l1] + H2*u[l2] + H3*u[l3];
00265 v[size2+i] = H3*u[l] - H2*u[l1] + H1*u[l2] - H0*u[l3];
00266 }
00267 for( int i = 0; i < size1; i++ ) {
00268 ii = i * size;
00269 wdata[ii+j] = v[i];
00270 }
00271 }
00272
00273 for( int i = 0; i < size1; i++ ) {
00274 ii = i * size;
00275 for( int j = 0; j < size1; j++ ) {
00276 u[j] = wdata[ii+j];
00277 }
00278 for( int j = 0; j < size2; j++ ) {
00279 l = 2 * j;
00280 l1 = (l+1)%size1;
00281 l2 = (l+2)%size1;
00282 l3 = (l+3)%size1;
00283 v[j] = H0*u[l] + H1*u[l1] + H2*u[l2] + H3*u[l3];
00284 v[size2+j] = H3*u[l] - H2*u[l1] + H1*u[l2] - H0*u[l3];
00285 }
00286 for( int j = 0; j < size1; j++ ) {
00287 wdata[ii+j] = v[j];
00288 }
00289 }
00290
00291 size1 = size2;
00292
00293 }
00294
00295 #ifdef STORM_FWT_CCHECKUSAGE
00296 if (size1 != ns) return 2;
00297 #endif
00298
00299 return 0;
00300
00301 #endif
00302 }
00303
00304
00305
00306 int pd4fwt2ds(
00307 const int& Sizef2d,
00308 const int& ns,
00309 const int& ks,
00310 const int& size,
00311 const double* sdata,
00312 double* wdata,
00313 double* ws
00314 ) {
00315
00316 #ifdef STORM_USEF
00317
00318 int iflag;
00319 dpd4fwt2dsws_( &Sizef2d, &ns, &ks, &size, sdata, wdata, ws, &iflag );
00320 return iflag;
00321
00322 #else
00323
00324 int ii,l;
00325 int l1,l2,l3;
00326 int size2;
00327 int size1 = size;
00328 double* u = ws;
00329 double* v = &ws[Sizef2d];
00330
00331 #ifdef STORM_FWT_CCHECKUSAGE
00332 if ( size1 < ns ) return 1;
00333 #endif
00334
00335 for( int i = 0; i < size1; i++ ) {
00336 ii = i * size;
00337 for( int j = 0; j < size1; j++ ) {
00338 wdata[ii+j] = sdata[ii+j];
00339 }
00340 }
00341
00342 if ( size1 == ns ) return 0;
00343
00344 while ( size1 > ns ) {
00345
00346 size2 = size1 / 2;
00347
00348 for( int j = 0; j < size1; j++ ) {
00349 for( int i = 0; i < size1; i++ ) {
00350 ii = i * size;
00351 u[i] = wdata[ii+j];
00352 }
00353 for( int i = 0; i < size2; i++ ) {
00354 l = 2 * i;
00355 l1 = (l+1)%size1;
00356 l2 = (l+2)%size1;
00357 l3 = (l+3)%size1;
00358 v[i] = H0*u[l] + H1*u[l1] + H2*u[l2] + H3*u[l3];
00359 v[size2+i] = H3*u[l] - H2*u[l1] + H1*u[l2] - H0*u[l3];
00360 }
00361 for( int i = 0; i < size1; i++ ) {
00362 ii = i * size;
00363 wdata[ii+j] = v[i];
00364 }
00365 }
00366
00367 for( int i = 0; i < size1; i++ ) {
00368 ii = i * size;
00369 for( int j = 0; j < size1; j++ ) {
00370 u[j] = wdata[ii+j];
00371 }
00372 for( int j = 0; j < size2; j++ ) {
00373 l = 2 * j;
00374 l1 = (l+1)%size1;
00375 l2 = (l+2)%size1;
00376 l3 = (l+3)%size1;
00377 v[j] = H0*u[l] + H1*u[l1] + H2*u[l2] + H3*u[l3];
00378 v[size2+j] = H3*u[l] - H2*u[l1] + H1*u[l2] - H0*u[l3];
00379 }
00380 for( int j = 0; j < size1; j++ ) {
00381 wdata[ii+j] = v[j];
00382 }
00383 }
00384
00385 size1 = size2;
00386
00387 }
00388
00389 #ifdef STORM_FWT_CCHECKUSAGE
00390 if (size1 != ns) return 2;
00391 #endif
00392
00393 return 0;
00394
00395 #endif
00396
00397 }
00398
00399
00400 #undef H0
00401 #undef H1
00402 #undef H2
00403 #undef H3
00404
