44 |
|
* |
45 |
|
* Created by Xiuquan Sun on 05/09/06. |
46 |
|
* @author Xiuquan Sun |
47 |
< |
* @version $Id: Hxy.cpp,v 1.2 2006-05-16 02:06:37 gezelter Exp $ |
47 |
> |
* @version $Id: Hxy.cpp,v 1.3 2006-05-16 20:38:23 gezelter Exp $ |
48 |
|
* |
49 |
|
*/ |
50 |
|
|
60 |
|
#include<string.h> |
61 |
|
#include<stdlib.h> |
62 |
|
#include<math.h> |
63 |
– |
#ifndef WITHOUT_FFTW |
64 |
– |
#include<fftw3.h> |
65 |
– |
#endif |
63 |
|
|
64 |
|
namespace oopse { |
65 |
|
|
87 |
|
} |
88 |
|
|
89 |
|
void Hxy::process() { |
90 |
< |
#ifndef WITHOUT_FFTW |
94 |
< |
|
90 |
> |
#if defined(HAVE_FFTW_H) || defined(HAVE_DFFTW_H) || defined(HAVE_FFTW3_H) |
91 |
|
DumpReader reader(info_, dumpFilename_); |
92 |
|
int nFrames = reader.getNFrames(); |
93 |
|
nProcessed_ = nFrames/step_; |
111 |
|
reader.readFrame(istep); |
112 |
|
currentSnapshot_ = info_->getSnapshotManager()->getCurrentSnapshot(); |
113 |
|
nMolecules = info_->getNGlobalMolecules(); |
114 |
< |
|
114 |
> |
|
115 |
|
Mat3x3d hmat = currentSnapshot_->getHmat(); |
120 |
– |
|
116 |
|
|
117 |
< |
fftw_complex *in, *out; |
117 |
> |
#ifdef HAVE_FFTW3_H |
118 |
|
fftw_plan p; |
119 |
+ |
#else |
120 |
+ |
fftwnd_plan p; |
121 |
+ |
#endif |
122 |
+ |
fftw_complex *in, *out; |
123 |
+ |
|
124 |
|
|
125 |
|
in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * (nBinsX_*nBinsY_)); |
126 |
|
out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) *(nBinsX_*nBinsY_)); |
127 |
< |
p = fftw_plan_dft_2d(nBinsX_, |
128 |
< |
nBinsY_, |
129 |
< |
in, out, |
130 |
< |
FFTW_FORWARD, |
131 |
< |
FFTW_ESTIMATE); |
127 |
> |
|
128 |
> |
#ifdef HAVE_FFTW3_H |
129 |
> |
p = fftw_plan_dft_2d(nBinsX_, nBinsY_, in, out, |
130 |
> |
FFTW_FORWARD, FFTW_ESTIMATE); |
131 |
> |
#else |
132 |
> |
p = fftw2d_create_plan(nBinsX_, nBinsY_, FFTW_FORWARD, FFTW_ESTIMATE); |
133 |
> |
#endif |
134 |
|
|
135 |
|
int i, j; |
136 |
|
|
267 |
|
c_im(in[newindex]) = 0.0; |
268 |
|
} |
269 |
|
} |
270 |
< |
|
271 |
< |
fftw_execute(p); |
270 |
> |
#ifdef HAVE_FFTW3_H |
271 |
> |
fftw_execute(p); |
272 |
> |
#else |
273 |
> |
fftwnd_one(p, in, out); |
274 |
> |
#endif |
275 |
|
|
276 |
|
for (i=0; i< nBinsX_; i++) { |
277 |
|
for(j=0; j< nBinsY_; j++) { |
279 |
|
mag[newindex] = pow(c_re(out[newindex]),2) + pow(c_im(out[newindex]),2); |
280 |
|
} |
281 |
|
} |
282 |
< |
|
282 |
> |
#ifdef HAVE_FFTW3_H |
283 |
|
fftw_destroy_plan(p); |
284 |
+ |
#else |
285 |
+ |
fftwnd_destroy_plan(p); |
286 |
+ |
#endif |
287 |
|
fftw_free(out); |
288 |
|
fftw_free(in); |
289 |
|
|