Source
xxxxxxxxxx
32
32
#include <casa/OS/HostInfo.h>
33
33
#include <synthesis/Utilities/FFT2D.h>
34
34
#include <lattices/Lattices/Lattice.h>
35
35
#ifdef _OPENMP
36
36
#include <omp.h>
37
37
#endif
38
38
using namespace casacore;
39
39
namespace casa { //# NAMESPACE CASA - BEGIN
40
40
41
41
// Utility function to try to give as much info as possible - CAS-12604
42
-
void raise_programmer_error(Long nx_p, Long ny_p, Long x, Long y, const char *file,
42
+
void throw_programmer_error(Long nx_p, Long ny_p, Long x, Long y, const char *file,
43
43
int line)
44
44
{
45
45
AipsError exc;
46
46
ostringstream msg;
47
47
msg << "Programmer error: called FFT2D with wrong size. In " << file << ":" << line
48
48
<< ", with nx_p: " << nx_p << ", ny_p: " << ny_p << ", x: " << x << ", y: " << y
49
49
<< "\n Stack trace: " << exc.getStackTrace();
50
50
exc.setMessage(msg.str());
51
51
throw(exc);
52
52
}
361
361
//Will need to seperate the plan from the execute if we want to run this in multiple threads
362
362
Int dim[2]={Int(y), Int(x)};
363
363
if(toFreq){
364
364
if(!planC2CD_forw_p){
365
365
planC2CD_forw_p=fftw_plan_dft(2, dim, reinterpret_cast<fftw_complex *>(out), reinterpret_cast<fftw_complex *>(out), FFTW_FORWARD, FFTW_ESTIMATE);
366
366
fftw_execute(planC2CD_forw_p);
367
367
nx_p=x;
368
368
ny_p=y;
369
369
}
370
370
else{
371
-
if(true or (nx_p != x) || (ny_p !=y)) {
372
-
raise_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
371
+
if((nx_p != x) || (ny_p !=y)) {
372
+
throw_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
373
373
}
374
374
fftw_execute_dft(planC2CD_forw_p, reinterpret_cast<fftw_complex *>(out), reinterpret_cast<fftw_complex *>(out));
375
375
}
376
376
//fft1_p.plan_c2c_forward(IPosition(2, x, y), out);
377
377
}
378
378
else{
379
379
if(!planC2CD_back_p){
380
380
planC2CD_back_p=fftw_plan_dft(2, dim, reinterpret_cast<fftw_complex *>(out), reinterpret_cast<fftw_complex *>(out), FFTW_BACKWARD, FFTW_ESTIMATE);
381
381
fftw_execute(planC2CD_back_p);
382
382
nx_p=x;
383
383
ny_p=y;
384
384
}
385
385
else{
386
-
if(true or (nx_p != x) || (ny_p !=y)) {
387
-
raise_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
386
+
if((nx_p != x) || (ny_p !=y)) {
387
+
throw_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
388
388
}
389
389
fftw_execute_dft(planC2CD_back_p, reinterpret_cast<fftw_complex *>(out), reinterpret_cast<fftw_complex *>(out));
390
390
}
391
391
// fft1_p.plan_c2c_backward(IPosition(2, x, y), out);
392
392
}
393
393
}
394
394
else{
395
395
throw(AipsError("Double precision FFT with FFTPack is not implemented"));
396
396
}
397
397
}
401
401
Int dim[2]={Int(y), Int(x)};
402
402
if(toFreq){
403
403
if(!planC2C_forw_p){
404
404
planC2C_forw_p=fftwf_plan_dft(2, dim, reinterpret_cast<fftwf_complex *>(out), reinterpret_cast<fftwf_complex *>(out), FFTW_FORWARD, FFTW_ESTIMATE);
405
405
fftwf_execute(planC2C_forw_p);
406
406
nx_p=x;
407
407
ny_p=y;
408
408
409
409
}
410
410
else{
411
-
if(true or (nx_p != x) || (ny_p !=y)) {
412
-
raise_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
411
+
if((nx_p != x) || (ny_p !=y)) {
412
+
throw_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
413
413
}
414
414
fftwf_execute_dft(planC2C_forw_p, reinterpret_cast<fftwf_complex *>(out), reinterpret_cast<fftwf_complex *>(out) );
415
415
}
416
416
//fft1_p.plan_c2c_forward(IPosition(2, x, y), out);
417
417
}
418
418
else{
419
419
if(!planC2C_back_p){
420
420
planC2C_back_p=fftwf_plan_dft(2, dim, reinterpret_cast<fftwf_complex *>(out), reinterpret_cast<fftwf_complex *>(out), FFTW_BACKWARD, FFTW_ESTIMATE);
421
421
fftwf_execute(planC2C_back_p);
422
422
nx_p=x;
423
423
ny_p=y;
424
424
}
425
425
else{
426
-
if(true or (nx_p != x) || (ny_p !=y)) {
427
-
raise_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
426
+
if((nx_p != x) || (ny_p !=y)) {
427
+
throw_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
428
428
}
429
429
fftwf_execute_dft(planC2C_back_p, reinterpret_cast<fftwf_complex *>(out), reinterpret_cast<fftwf_complex *>(out) );
430
430
}
431
431
// fft1_p.plan_c2c_backward(IPosition(2, x, y), out);
432
432
}
433
433
434
434
435
435
}
436
436
else{
437
437
Int ier;
459
459
if(!planR2C_p){
460
460
planR2C_p=fftwf_plan_dft_r2c(2, dim, in, reinterpret_cast<fftwf_complex *>(out), FFTW_ESTIMATE);
461
461
462
462
//fft1_p.plan_c2c_forward(IPosition(2, x, y), out);
463
463
464
464
fftwf_execute(planR2C_p);
465
465
nx_p=x;
466
466
ny_p=y;
467
467
}
468
468
else{
469
-
if(true or (nx_p != x) || (ny_p !=y)) {
470
-
raise_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
469
+
if((nx_p != x) || (ny_p !=y)) {
470
+
throw_programmer_error(nx_p, ny_p, x, y, __FILE__, __LINE__);
471
471
}
472
472
fftwf_execute_dft_r2c(planR2C_p, in, reinterpret_cast<fftwf_complex *>(out));
473
473
}
474
474
}
475
475
else{
476
476
/*
477
477
Int ier;
478
478
Int x1=Int(x);
479
479
Int y1=Int(y);
480
480
if(wsave_p.size()==0){