Mercurial > hg > octave-lyh

annotate liboctave/oct-fftw.cc @ 11188:4cb1522e4d0f

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Use function handle as input to cellfun, rather than quoted function name or anonymous function wrapper.
author Rik <octave@nomad.inbox5.com>
date Wed, 03 Nov 2010 17:20:56 -0700
parents 12884915a8e4
children 331fcc41ca23
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
1 /*
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
2
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
3 Copyright (C) 2001, 2002, 2004, 2005, 2006, 2007, 2008, 2009 John W. Eaton
7017
a1dbe9d80eee [project @ 2007-10-12 21:27:11 by jwe]
jwe
parents: 7016
diff changeset
4
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
5 This file is part of Octave.
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
6
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
7 Octave is free software; you can redistribute it and/or modify it
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
8 under the terms of the GNU General Public License as published by the
7016
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 6482
diff changeset
9 Free Software Foundation; either version 3 of the License, or (at your
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 6482
diff changeset
10 option) any later version.
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
11
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
12 Octave is distributed in the hope that it will be useful, but WITHOUT
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
15 for more details.
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
16
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
17 You should have received a copy of the GNU General Public License
7016
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 6482
diff changeset
18 along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 6482
diff changeset
19 <http://www.gnu.org/licenses/>.
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
20
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
21 */
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
22
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
23 #ifdef HAVE_CONFIG_H
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
24 #include <config.h>
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
25 #endif
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
26
9523
0ce82753dd72 more configure changes for libraries
John W. Eaton <jwe@octave.org>
parents: 9516
diff changeset
27 #if defined (HAVE_FFTW)
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
28
4775
88b638195bd1 [project @ 2004-02-16 19:57:06 by jwe]
jwe
parents: 4774
diff changeset
29 #include <iostream>
88b638195bd1 [project @ 2004-02-16 19:57:06 by jwe]
jwe
parents: 4774
diff changeset
30 #include <vector>
88b638195bd1 [project @ 2004-02-16 19:57:06 by jwe]
jwe
parents: 4774
diff changeset
31
4786
fc316bde0053 [project @ 2004-02-18 12:52:20 by jwe]
jwe
parents: 4785
diff changeset
32 #include "lo-error.h"
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
33 #include "oct-fftw.h"
4786
fc316bde0053 [project @ 2004-02-18 12:52:20 by jwe]
jwe
parents: 4785
diff changeset
34 #include "quit.h"
8377
25bc2d31e1bf improve OCTAVE_LOCAL_BUFFER
Jaroslav Hajek <highegg@gmail.com>
parents: 7789
diff changeset
35 #include "oct-locbuf.h"
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
36
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
37 octave_fftw_planner *octave_fftw_planner::instance = 0;
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
38
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
39 // Helper class to create and cache fftw plans for both 1d and
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
40 // 2d. This implementation defaults to using FFTW_ESTIMATE to create
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
41 // the plans, which in theory is suboptimal, but provides quit
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
42 // reasonable performance.
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
43
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
44 // Also note that if FFTW_ESTIMATE is not used the planner in FFTW3
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
45 // destroys the input and output arrays. We must therefore create a
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
46 // temporary input array with the same size and 16-byte alignment as
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
47 // the original array and use that for the planner. Note that we also
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
48 // use any wisdom that is available, either in a FFTW3 system wide file
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
49 // or as supplied by the user.
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
50
5775
ace8d8d26933 [project @ 2006-04-24 19:13:06 by jwe]
jwe
parents: 5307
diff changeset
51 // FIXME -- if we can ensure 16 byte alignment in Array<T>
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
52 // (<T> *data) the FFTW3 can use SIMD instructions for further
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
53 // acceleration.
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
54
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
55 // Note that it is profitable to store the FFTW3 plans, for small
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
56 // ffts.
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
57
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
58 octave_fftw_planner::octave_fftw_planner (void)
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
59 {
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
60 meth = ESTIMATE;
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
61
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
62 plan[0] = plan[1] = 0;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
63 d[0] = d[1] = s[0] = s[1] = r[0] = r[1] = h[0] = h[1] = 0;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
64 simd_align[0] = simd_align[1] = false;
5044
50140aa23b32 [project @ 2004-10-11 16:31:51 by jwe]
jwe
parents: 4809
diff changeset
65 inplace[0] = inplace[1] = false;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
66 n[0] = n[1] = dim_vector ();
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
67
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
68 rplan = 0;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
69 rd = rs = rr = rh = 0;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
70 rsimd_align = false;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
71 rn = dim_vector ();
5044
50140aa23b32 [project @ 2004-10-11 16:31:51 by jwe]
jwe
parents: 4809
diff changeset
72
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
73 // If we have a system wide wisdom file, import it.
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
74 fftw_import_system_wisdom ();
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
75 }
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
76
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
77 bool
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
78 octave_fftw_planner::instance_ok (void)
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
79 {
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
80 bool retval = true;
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
81
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
82 if (! instance)
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
83 instance = new octave_fftw_planner ();
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
84
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
85 if (! instance)
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
86 {
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
87 (*current_liboctave_error_handler)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
88 ("unable to create octave_fftw_planner object!");
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
89
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
90 retval = false;
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
91 }
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
92
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
93 return retval;
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
94 }
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
95
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
96 #define CHECK_SIMD_ALIGNMENT(x) \
6482
7e958a1532c6 [project @ 2007-04-04 02:37:07 by jwe]
jwe
parents: 6228
diff changeset
97 (((reinterpret_cast<ptrdiff_t> (x)) & 0xF) == 0)
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
98
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
99 fftw_plan
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
100 octave_fftw_planner::do_create_plan (int dir, const int rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
101 const dim_vector dims,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
102 octave_idx_type howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
103 octave_idx_type stride,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
104 octave_idx_type dist,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
105 const Complex *in, Complex *out)
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
106 {
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
107 int which = (dir == FFTW_FORWARD) ? 0 : 1;
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
108 fftw_plan *cur_plan_p = &plan[which];
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
109 bool create_new_plan = false;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
110 bool ioalign = CHECK_SIMD_ALIGNMENT (in) && CHECK_SIMD_ALIGNMENT (out);
5044
50140aa23b32 [project @ 2004-10-11 16:31:51 by jwe]
jwe
parents: 4809
diff changeset
111 bool ioinplace = (in == out);
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
112
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
113 // Don't create a new plan if we have a non SIMD plan already but
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
114 // can do SIMD. This prevents endlessly recreating plans if we
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
115 // change the alignment.
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
116
4783
55975a3073be [project @ 2004-02-17 17:16:37 by jwe]
jwe
parents: 4775
diff changeset
117 if (plan[which] == 0 || d[which] != dist || s[which] != stride
5044
50140aa23b32 [project @ 2004-10-11 16:31:51 by jwe]
jwe
parents: 4809
diff changeset
118 || r[which] != rank || h[which] != howmany
50140aa23b32 [project @ 2004-10-11 16:31:51 by jwe]
jwe
parents: 4809
diff changeset
119 || ioinplace != inplace[which]
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
120 || ((ioalign != simd_align[which]) ? !ioalign : false))
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
121 create_new_plan = true;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
122 else
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
123 {
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
124 // We still might not have the same shape of array.
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
125
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
126 for (int i = 0; i < rank; i++)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
127 if (dims(i) != n[which](i))
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
128 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
129 create_new_plan = true;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
130 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
131 }
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
132 }
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
133
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
134 if (create_new_plan)
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
135 {
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
136 d[which] = dist;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
137 s[which] = stride;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
138 r[which] = rank;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
139 h[which] = howmany;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
140 simd_align[which] = ioalign;
5044
50140aa23b32 [project @ 2004-10-11 16:31:51 by jwe]
jwe
parents: 4809
diff changeset
141 inplace[which] = ioinplace;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
142 n[which] = dims;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
143
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
144 // Note reversal of dimensions for column major storage in FFTW.
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
145 octave_idx_type nn = 1;
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
146 OCTAVE_LOCAL_BUFFER (int, tmp, rank);
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
147
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
148 for (int i = 0, j = rank-1; i < rank; i++, j--)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
149 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
150 tmp[i] = dims(j);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
151 nn *= dims(j);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
152 }
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
153
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
154 int plan_flags = 0;
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
155 bool plan_destroys_in = true;
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
156
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
157 switch (meth)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
158 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
159 case UNKNOWN:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
160 case ESTIMATE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
161 plan_flags |= FFTW_ESTIMATE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
162 plan_destroys_in = false;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
163 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
164 case MEASURE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
165 plan_flags |= FFTW_MEASURE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
166 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
167 case PATIENT:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
168 plan_flags |= FFTW_PATIENT;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
169 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
170 case EXHAUSTIVE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
171 plan_flags |= FFTW_EXHAUSTIVE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
172 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
173 case HYBRID:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
174 if (nn < 8193)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
175 plan_flags |= FFTW_MEASURE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
176 else
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
177 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
178 plan_flags |= FFTW_ESTIMATE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
179 plan_destroys_in = false;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
180 }
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
181 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
182 }
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
183
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
184 if (ioalign)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
185 plan_flags &= ~FFTW_UNALIGNED;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
186 else
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
187 plan_flags |= FFTW_UNALIGNED;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
188
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
189 if (*cur_plan_p)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
190 fftw_destroy_plan (*cur_plan_p);
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
191
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
192 if (plan_destroys_in)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
193 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
194 // Create matrix with the same size and 16-byte alignment as input
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
195 OCTAVE_LOCAL_BUFFER (Complex, itmp, nn * howmany + 32);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
196 itmp = reinterpret_cast<Complex *>
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
197 (((reinterpret_cast<ptrdiff_t>(itmp) + 15) & ~ 0xF) +
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
198 ((reinterpret_cast<ptrdiff_t> (in)) & 0xF));
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
199
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
200 *cur_plan_p =
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
201 fftw_plan_many_dft (rank, tmp, howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
202 reinterpret_cast<fftw_complex *> (itmp),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
203 0, stride, dist, reinterpret_cast<fftw_complex *> (out),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
204 0, stride, dist, dir, plan_flags);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
205 }
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
206 else
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
207 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
208 *cur_plan_p =
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
209 fftw_plan_many_dft (rank, tmp, howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
210 reinterpret_cast<fftw_complex *> (const_cast<Complex *> (in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
211 0, stride, dist, reinterpret_cast<fftw_complex *> (out),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
212 0, stride, dist, dir, plan_flags);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
213 }
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
214
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
215 if (*cur_plan_p == 0)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
216 (*current_liboctave_error_handler) ("Error creating fftw plan");
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
217 }
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
218
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
219 return *cur_plan_p;
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
220 }
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
221
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
222 fftw_plan
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
223 octave_fftw_planner::do_create_plan (const int rank, const dim_vector dims,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
224 octave_idx_type howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
225 octave_idx_type stride,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
226 octave_idx_type dist,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
227 const double *in, Complex *out)
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
228 {
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
229 fftw_plan *cur_plan_p = &rplan;
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
230 bool create_new_plan = false;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
231 bool ioalign = CHECK_SIMD_ALIGNMENT (in) && CHECK_SIMD_ALIGNMENT (out);
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
232
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
233 // Don't create a new plan if we have a non SIMD plan already but
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
234 // can do SIMD. This prevents endlessly recreating plans if we
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
235 // change the alignment.
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
236
4783
55975a3073be [project @ 2004-02-17 17:16:37 by jwe]
jwe
parents: 4775
diff changeset
237 if (rplan == 0 || rd != dist || rs != stride || rr != rank
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
238 || rh != howmany || ((ioalign != rsimd_align) ? !ioalign : false))
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
239 create_new_plan = true;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
240 else
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
241 {
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
242 // We still might not have the same shape of array.
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
243
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
244 for (int i = 0; i < rank; i++)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
245 if (dims(i) != rn(i))
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
246 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
247 create_new_plan = true;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
248 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
249 }
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
250 }
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
251
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
252 if (create_new_plan)
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
253 {
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
254 rd = dist;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
255 rs = stride;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
256 rr = rank;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
257 rh = howmany;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
258 rsimd_align = ioalign;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
259 rn = dims;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
260
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
261 // Note reversal of dimensions for column major storage in FFTW.
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
262 octave_idx_type nn = 1;
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
263 OCTAVE_LOCAL_BUFFER (int, tmp, rank);
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
264
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
265 for (int i = 0, j = rank-1; i < rank; i++, j--)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
266 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
267 tmp[i] = dims(j);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
268 nn *= dims(j);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
269 }
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
270
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
271 int plan_flags = 0;
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
272 bool plan_destroys_in = true;
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
273
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
274 switch (meth)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
275 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
276 case UNKNOWN:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
277 case ESTIMATE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
278 plan_flags |= FFTW_ESTIMATE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
279 plan_destroys_in = false;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
280 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
281 case MEASURE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
282 plan_flags |= FFTW_MEASURE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
283 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
284 case PATIENT:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
285 plan_flags |= FFTW_PATIENT;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
286 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
287 case EXHAUSTIVE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
288 plan_flags |= FFTW_EXHAUSTIVE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
289 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
290 case HYBRID:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
291 if (nn < 8193)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
292 plan_flags |= FFTW_MEASURE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
293 else
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
294 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
295 plan_flags |= FFTW_ESTIMATE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
296 plan_destroys_in = false;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
297 }
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
298 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
299 }
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
300
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
301 if (ioalign)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
302 plan_flags &= ~FFTW_UNALIGNED;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
303 else
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
304 plan_flags |= FFTW_UNALIGNED;
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
305
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
306 if (*cur_plan_p)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
307 fftw_destroy_plan (*cur_plan_p);
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
308
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
309 if (plan_destroys_in)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
310 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
311 // Create matrix with the same size and 16-byte alignment as input
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
312 OCTAVE_LOCAL_BUFFER (double, itmp, nn + 32);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
313 itmp = reinterpret_cast<double *>
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
314 (((reinterpret_cast<ptrdiff_t>(itmp) + 15) & ~ 0xF) +
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
315 ((reinterpret_cast<ptrdiff_t> (in)) & 0xF));
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
316
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
317 *cur_plan_p =
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
318 fftw_plan_many_dft_r2c (rank, tmp, howmany, itmp,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
319 0, stride, dist, reinterpret_cast<fftw_complex *> (out),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
320 0, stride, dist, plan_flags);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
321 }
6228
aa5df9ba98d5 [project @ 2007-01-05 22:49:03 by dbateman]
dbateman
parents: 5775
diff changeset
322 else
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
323 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
324 *cur_plan_p =
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
325 fftw_plan_many_dft_r2c (rank, tmp, howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
326 (const_cast<double *> (in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
327 0, stride, dist, reinterpret_cast<fftw_complex *> (out),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
328 0, stride, dist, plan_flags);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
329 }
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
330
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
331 if (*cur_plan_p == 0)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
332 (*current_liboctave_error_handler) ("Error creating fftw plan");
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
333 }
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
334
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
335 return *cur_plan_p;
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
336 }
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
337
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
338 octave_fftw_planner::FftwMethod
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
339 octave_fftw_planner::do_method (void)
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
340 {
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
341 return meth;
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
342 }
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
343
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
344 octave_fftw_planner::FftwMethod
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
345 octave_fftw_planner::do_method (FftwMethod _meth)
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
346 {
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
347 FftwMethod ret = meth;
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
348 if (_meth == ESTIMATE || _meth == MEASURE
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
349 || _meth == PATIENT || _meth == EXHAUSTIVE
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
350 || _meth == HYBRID)
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
351 {
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
352 if (meth != _meth)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
353 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
354 meth = _meth;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
355 if (rplan)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
356 fftw_destroy_plan (rplan);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
357 if (plan[0])
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
358 fftw_destroy_plan (plan[0]);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
359 if (plan[1])
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
360 fftw_destroy_plan (plan[1]);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
361 rplan = plan[0] = plan[1] = 0;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
362 }
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
363 }
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
364 else
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
365 ret = UNKNOWN;
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
366 return ret;
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
367 }
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
368
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
369 octave_float_fftw_planner *octave_float_fftw_planner::instance = 0;
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
370
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
371 octave_float_fftw_planner::octave_float_fftw_planner (void)
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
372 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
373 meth = ESTIMATE;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
374
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
375 plan[0] = plan[1] = 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
376 d[0] = d[1] = s[0] = s[1] = r[0] = r[1] = h[0] = h[1] = 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
377 simd_align[0] = simd_align[1] = false;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
378 inplace[0] = inplace[1] = false;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
379 n[0] = n[1] = dim_vector ();
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
380
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
381 rplan = 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
382 rd = rs = rr = rh = 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
383 rsimd_align = false;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
384 rn = dim_vector ();
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
385
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
386 // If we have a system wide wisdom file, import it.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
387 fftwf_import_system_wisdom ();
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
388 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
389
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
390 bool
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
391 octave_float_fftw_planner::instance_ok (void)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
392 {
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
393 bool retval = true;
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
394
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
395 if (! instance)
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
396 instance = new octave_float_fftw_planner ();
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
397
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
398 if (! instance)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
399 {
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
400 (*current_liboctave_error_handler)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
401 ("unable to create octave_fftw_planner object!");
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
402
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
403 retval = false;
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
404 }
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
405
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
406 return retval;
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
407 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
408
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
409 fftwf_plan
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
410 octave_float_fftw_planner::do_create_plan (int dir, const int rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
411 const dim_vector dims,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
412 octave_idx_type howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
413 octave_idx_type stride,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
414 octave_idx_type dist,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
415 const FloatComplex *in,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
416 FloatComplex *out)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
417 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
418 int which = (dir == FFTW_FORWARD) ? 0 : 1;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
419 fftwf_plan *cur_plan_p = &plan[which];
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
420 bool create_new_plan = false;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
421 bool ioalign = CHECK_SIMD_ALIGNMENT (in) && CHECK_SIMD_ALIGNMENT (out);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
422 bool ioinplace = (in == out);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
423
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
424 // Don't create a new plan if we have a non SIMD plan already but
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
425 // can do SIMD. This prevents endlessly recreating plans if we
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
426 // change the alignment.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
427
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
428 if (plan[which] == 0 || d[which] != dist || s[which] != stride
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
429 || r[which] != rank || h[which] != howmany
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
430 || ioinplace != inplace[which]
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
431 || ((ioalign != simd_align[which]) ? !ioalign : false))
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
432 create_new_plan = true;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
433 else
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
434 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
435 // We still might not have the same shape of array.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
436
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
437 for (int i = 0; i < rank; i++)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
438 if (dims(i) != n[which](i))
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
439 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
440 create_new_plan = true;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
441 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
442 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
443 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
444
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
445 if (create_new_plan)
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
446 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
447 d[which] = dist;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
448 s[which] = stride;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
449 r[which] = rank;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
450 h[which] = howmany;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
451 simd_align[which] = ioalign;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
452 inplace[which] = ioinplace;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
453 n[which] = dims;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
454
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
455 // Note reversal of dimensions for column major storage in FFTW.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
456 octave_idx_type nn = 1;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
457 OCTAVE_LOCAL_BUFFER (int, tmp, rank);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
458
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
459 for (int i = 0, j = rank-1; i < rank; i++, j--)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
460 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
461 tmp[i] = dims(j);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
462 nn *= dims(j);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
463 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
464
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
465 int plan_flags = 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
466 bool plan_destroys_in = true;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
467
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
468 switch (meth)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
469 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
470 case UNKNOWN:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
471 case ESTIMATE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
472 plan_flags |= FFTW_ESTIMATE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
473 plan_destroys_in = false;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
474 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
475 case MEASURE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
476 plan_flags |= FFTW_MEASURE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
477 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
478 case PATIENT:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
479 plan_flags |= FFTW_PATIENT;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
480 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
481 case EXHAUSTIVE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
482 plan_flags |= FFTW_EXHAUSTIVE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
483 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
484 case HYBRID:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
485 if (nn < 8193)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
486 plan_flags |= FFTW_MEASURE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
487 else
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
488 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
489 plan_flags |= FFTW_ESTIMATE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
490 plan_destroys_in = false;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
491 }
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
492 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
493 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
494
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
495 if (ioalign)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
496 plan_flags &= ~FFTW_UNALIGNED;
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
497 else
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
498 plan_flags |= FFTW_UNALIGNED;
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
499
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
500 if (*cur_plan_p)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
501 fftwf_destroy_plan (*cur_plan_p);
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
502
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
503 if (plan_destroys_in)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
504 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
505 // Create matrix with the same size and 16-byte alignment as input
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
506 OCTAVE_LOCAL_BUFFER (FloatComplex, itmp, nn * howmany + 32);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
507 itmp = reinterpret_cast<FloatComplex *>
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
508 (((reinterpret_cast<ptrdiff_t>(itmp) + 15) & ~ 0xF) +
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
509 ((reinterpret_cast<ptrdiff_t> (in)) & 0xF));
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
510
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
511 *cur_plan_p =
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
512 fftwf_plan_many_dft (rank, tmp, howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
513 reinterpret_cast<fftwf_complex *> (itmp),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
514 0, stride, dist, reinterpret_cast<fftwf_complex *> (out),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
515 0, stride, dist, dir, plan_flags);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
516 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
517 else
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
518 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
519 *cur_plan_p =
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
520 fftwf_plan_many_dft (rank, tmp, howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
521 reinterpret_cast<fftwf_complex *> (const_cast<FloatComplex *> (in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
522 0, stride, dist, reinterpret_cast<fftwf_complex *> (out),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
523 0, stride, dist, dir, plan_flags);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
524 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
525
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
526 if (*cur_plan_p == 0)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
527 (*current_liboctave_error_handler) ("Error creating fftw plan");
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
528 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
529
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
530 return *cur_plan_p;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
531 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
532
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
533 fftwf_plan
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
534 octave_float_fftw_planner::do_create_plan (const int rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
535 const dim_vector dims,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
536 octave_idx_type howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
537 octave_idx_type stride,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
538 octave_idx_type dist,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
539 const float *in, FloatComplex *out)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
540 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
541 fftwf_plan *cur_plan_p = &rplan;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
542 bool create_new_plan = false;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
543 bool ioalign = CHECK_SIMD_ALIGNMENT (in) && CHECK_SIMD_ALIGNMENT (out);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
544
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
545 // Don't create a new plan if we have a non SIMD plan already but
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
546 // can do SIMD. This prevents endlessly recreating plans if we
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
547 // change the alignment.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
548
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
549 if (rplan == 0 || rd != dist || rs != stride || rr != rank
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
550 || rh != howmany || ((ioalign != rsimd_align) ? !ioalign : false))
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
551 create_new_plan = true;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
552 else
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
553 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
554 // We still might not have the same shape of array.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
555
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
556 for (int i = 0; i < rank; i++)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
557 if (dims(i) != rn(i))
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
558 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
559 create_new_plan = true;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
560 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
561 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
562 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
563
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
564 if (create_new_plan)
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
565 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
566 rd = dist;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
567 rs = stride;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
568 rr = rank;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
569 rh = howmany;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
570 rsimd_align = ioalign;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
571 rn = dims;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
572
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
573 // Note reversal of dimensions for column major storage in FFTW.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
574 octave_idx_type nn = 1;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
575 OCTAVE_LOCAL_BUFFER (int, tmp, rank);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
576
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
577 for (int i = 0, j = rank-1; i < rank; i++, j--)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
578 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
579 tmp[i] = dims(j);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
580 nn *= dims(j);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
581 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
582
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
583 int plan_flags = 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
584 bool plan_destroys_in = true;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
585
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
586 switch (meth)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
587 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
588 case UNKNOWN:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
589 case ESTIMATE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
590 plan_flags |= FFTW_ESTIMATE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
591 plan_destroys_in = false;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
592 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
593 case MEASURE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
594 plan_flags |= FFTW_MEASURE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
595 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
596 case PATIENT:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
597 plan_flags |= FFTW_PATIENT;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
598 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
599 case EXHAUSTIVE:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
600 plan_flags |= FFTW_EXHAUSTIVE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
601 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
602 case HYBRID:
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
603 if (nn < 8193)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
604 plan_flags |= FFTW_MEASURE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
605 else
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
606 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
607 plan_flags |= FFTW_ESTIMATE;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
608 plan_destroys_in = false;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
609 }
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
610 break;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
611 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
612
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
613 if (ioalign)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
614 plan_flags &= ~FFTW_UNALIGNED;
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
615 else
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
616 plan_flags |= FFTW_UNALIGNED;
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
617
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
618 if (*cur_plan_p)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
619 fftwf_destroy_plan (*cur_plan_p);
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
620
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
621 if (plan_destroys_in)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
622 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
623 // Create matrix with the same size and 16-byte alignment as input
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
624 OCTAVE_LOCAL_BUFFER (float, itmp, nn + 32);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
625 itmp = reinterpret_cast<float *>
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
626 (((reinterpret_cast<ptrdiff_t>(itmp) + 15) & ~ 0xF) +
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
627 ((reinterpret_cast<ptrdiff_t> (in)) & 0xF));
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
628
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
629 *cur_plan_p =
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
630 fftwf_plan_many_dft_r2c (rank, tmp, howmany, itmp,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
631 0, stride, dist, reinterpret_cast<fftwf_complex *> (out),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
632 0, stride, dist, plan_flags);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
633 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
634 else
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
635 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
636 *cur_plan_p =
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
637 fftwf_plan_many_dft_r2c (rank, tmp, howmany,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
638 (const_cast<float *> (in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
639 0, stride, dist, reinterpret_cast<fftwf_complex *> (out),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
640 0, stride, dist, plan_flags);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
641 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
642
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
643 if (*cur_plan_p == 0)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
644 (*current_liboctave_error_handler) ("Error creating fftw plan");
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
645 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
646
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
647 return *cur_plan_p;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
648 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
649
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
650 octave_float_fftw_planner::FftwMethod
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
651 octave_float_fftw_planner::do_method (void)
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
652 {
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
653 return meth;
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
654 }
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
655
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
656 octave_float_fftw_planner::FftwMethod
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
657 octave_float_fftw_planner::do_method (FftwMethod _meth)
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
658 {
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
659 FftwMethod ret = meth;
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
660 if (_meth == ESTIMATE || _meth == MEASURE
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
661 || _meth == PATIENT || _meth == EXHAUSTIVE
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
662 || _meth == HYBRID)
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
663 {
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
664 if (meth != _meth)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
665 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
666 meth = _meth;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
667 if (rplan)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
668 fftwf_destroy_plan (rplan);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
669 if (plan[0])
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
670 fftwf_destroy_plan (plan[0]);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
671 if (plan[1])
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
672 fftwf_destroy_plan (plan[1]);
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
673 rplan = plan[0] = plan[1] = 0;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
674 }
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
675 }
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
676 else
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
677 ret = UNKNOWN;
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
678 return ret;
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
679 }
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
680
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
681 template <class T>
4775
88b638195bd1 [project @ 2004-02-16 19:57:06 by jwe]
jwe
parents: 4774
diff changeset
682 static inline void
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
683 convert_packcomplex_1d (T *out, size_t nr, size_t nc,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
684 octave_idx_type stride, octave_idx_type dist)
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
685 {
10142
829e69ec3110 make OCTAVE_QUIT a function
Jaroslav Hajek <highegg@gmail.com>
parents: 9523
diff changeset
686 octave_quit ();
4785
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
687
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
688 // Fill in the missing data.
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
689
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
690 for (size_t i = 0; i < nr; i++)
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
691 for (size_t j = nc/2+1; j < nc; j++)
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
692 out[j*stride + i*dist] = conj(out[(nc - j)*stride + i*dist]);
4785
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
693
10142
829e69ec3110 make OCTAVE_QUIT a function
Jaroslav Hajek <highegg@gmail.com>
parents: 9523
diff changeset
694 octave_quit ();
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
695 }
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
696
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
697 template <class T>
4775
88b638195bd1 [project @ 2004-02-16 19:57:06 by jwe]
jwe
parents: 4774
diff changeset
698 static inline void
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
699 convert_packcomplex_Nd (T *out, const dim_vector &dv)
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
700 {
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
701 size_t nc = dv(0);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
702 size_t nr = dv(1);
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
703 size_t np = (dv.length () > 2 ? dv.numel () / nc / nr : 1);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
704 size_t nrp = nr * np;
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
705 T *ptr1, *ptr2;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
706
10142
829e69ec3110 make OCTAVE_QUIT a function
Jaroslav Hajek <highegg@gmail.com>
parents: 9523
diff changeset
707 octave_quit ();
4785
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
708
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
709 // Create space for the missing elements.
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
710
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
711 for (size_t i = 0; i < nrp; i++)
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
712 {
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
713 ptr1 = out + i * (nc/2 + 1) + nrp*((nc-1)/2);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
714 ptr2 = out + i * nc;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
715 for (size_t j = 0; j < nc/2+1; j++)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
716 *ptr2++ = *ptr1++;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
717 }
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
718
10142
829e69ec3110 make OCTAVE_QUIT a function
Jaroslav Hajek <highegg@gmail.com>
parents: 9523
diff changeset
719 octave_quit ();
4785
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
720
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
721 // Fill in the missing data for the rank = 2 case directly for speed.
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
722
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
723 for (size_t i = 0; i < np; i++)
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
724 {
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
725 for (size_t j = 1; j < nr; j++)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
726 for (size_t k = nc/2+1; k < nc; k++)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
727 out[k + (j + i*nr)*nc] = conj(out[nc - k + ((i+1)*nr - j)*nc]);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
728
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
729 for (size_t j = nc/2+1; j < nc; j++)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
730 out[j + i*nr*nc] = conj(out[(i*nr+1)*nc - j]);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
731 }
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
732
10142
829e69ec3110 make OCTAVE_QUIT a function
Jaroslav Hajek <highegg@gmail.com>
parents: 9523
diff changeset
733 octave_quit ();
4785
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
734
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
735 // Now do the permutations needed for rank > 2 cases.
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
736
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
737 size_t jstart = dv(0) * dv(1);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
738 size_t kstep = dv(0);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
739 size_t nel = dv.numel ();
4785
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
740
4808
a9ec0ce18568 [project @ 2004-03-02 17:36:28 by jwe]
jwe
parents: 4788
diff changeset
741 for (int inner = 2; inner < dv.length (); inner++)
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
742 {
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
743 size_t jmax = jstart * dv(inner);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
744 for (size_t i = 0; i < nel; i+=jmax)
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
745 for (size_t j = jstart, jj = jmax-jstart; j < jj;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
746 j+=jstart, jj-=jstart)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
747 for (size_t k = 0; k < jstart; k+= kstep)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
748 for (size_t l = nc/2+1; l < nc; l++)
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
749 {
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
750 T tmp = out[i+ j + k + l];
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
751 out[i + j + k + l] = out[i + jj + k + l];
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
752 out[i + jj + k + l] = tmp;
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
753 }
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
754 jstart = jmax;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
755 }
4785
d3018a33c584 [project @ 2004-02-17 21:58:17 by jwe]
jwe
parents: 4784
diff changeset
756
10142
829e69ec3110 make OCTAVE_QUIT a function
Jaroslav Hajek <highegg@gmail.com>
parents: 9523
diff changeset
757 octave_quit ();
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
758 }
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
759
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
760 int
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
761 octave_fftw::fft (const double *in, Complex *out, size_t npts,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
762 size_t nsamples, octave_idx_type stride, octave_idx_type dist)
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
763 {
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
764 dist = (dist < 0 ? npts : dist);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
765
10350
12884915a8e4 merge MArray classes & improve Array interface
Jaroslav Hajek <highegg@gmail.com>
parents: 10314
diff changeset
766 dim_vector dv (npts, 1);
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
767 fftw_plan plan = octave_fftw_planner::create_plan (1, dv, nsamples,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
768 stride, dist, in, out);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
769
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
770 fftw_execute_dft_r2c (plan, (const_cast<double *>(in)),
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
771 reinterpret_cast<fftw_complex *> (out));
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
772
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
773 // Need to create other half of the transform.
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
774
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
775 convert_packcomplex_1d (out, nsamples, npts, stride, dist);
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
776
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
777 return 0;
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
778 }
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
779
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
780 int
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
781 octave_fftw::fft (const Complex *in, Complex *out, size_t npts,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
782 size_t nsamples, octave_idx_type stride, octave_idx_type dist)
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
783 {
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
784 dist = (dist < 0 ? npts : dist);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
785
10350
12884915a8e4 merge MArray classes & improve Array interface
Jaroslav Hajek <highegg@gmail.com>
parents: 10314
diff changeset
786 dim_vector dv (npts, 1);
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
787 fftw_plan plan = octave_fftw_planner::create_plan (FFTW_FORWARD, 1, dv,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
788 nsamples, stride,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
789 dist, in, out);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
790
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
791 fftw_execute_dft (plan,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
792 reinterpret_cast<fftw_complex *> (const_cast<Complex *>(in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
793 reinterpret_cast<fftw_complex *> (out));
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
794
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
795 return 0;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
796 }
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
797
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
798 int
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
799 octave_fftw::ifft (const Complex *in, Complex *out, size_t npts,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
800 size_t nsamples, octave_idx_type stride, octave_idx_type dist)
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
801 {
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
802 dist = (dist < 0 ? npts : dist);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
803
10350
12884915a8e4 merge MArray classes & improve Array interface
Jaroslav Hajek <highegg@gmail.com>
parents: 10314
diff changeset
804 dim_vector dv (npts, 1);
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
805 fftw_plan plan = octave_fftw_planner::create_plan (FFTW_BACKWARD, 1, dv,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
806 nsamples, stride,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
807 dist, in, out);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
808
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
809 fftw_execute_dft (plan,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
810 reinterpret_cast<fftw_complex *> (const_cast<Complex *>(in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
811 reinterpret_cast<fftw_complex *> (out));
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
812
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
813 const Complex scale = npts;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
814 for (size_t j = 0; j < nsamples; j++)
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
815 for (size_t i = 0; i < npts; i++)
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
816 out[i*stride + j*dist] /= scale;
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
817
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
818 return 0;
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
819 }
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
820
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
821 int
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
822 octave_fftw::fftNd (const double *in, Complex *out, const int rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
823 const dim_vector &dv)
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
824 {
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5044
diff changeset
825 octave_idx_type dist = 1;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
826 for (int i = 0; i < rank; i++)
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
827 dist *= dv(i);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
828
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
829 // Fool with the position of the start of the output matrix, so that
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
830 // creating other half of the matrix won't cause cache problems.
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
831
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5044
diff changeset
832 octave_idx_type offset = (dv.numel () / dv(0)) * ((dv(0) - 1) / 2);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
833
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
834 fftw_plan plan = octave_fftw_planner::create_plan (rank, dv, 1, 1, dist,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
835 in, out + offset);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
836
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
837 fftw_execute_dft_r2c (plan, (const_cast<double *>(in)),
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
838 reinterpret_cast<fftw_complex *> (out+ offset));
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
839
4809
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
840 // Need to create other half of the transform.
b60be7678bdc [project @ 2004-03-02 17:40:08 by jwe]
jwe
parents: 4808
diff changeset
841
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
842 convert_packcomplex_Nd (out, dv);
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
843
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
844 return 0;
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
845 }
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
846
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
847 int
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
848 octave_fftw::fftNd (const Complex *in, Complex *out, const int rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
849 const dim_vector &dv)
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
850 {
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5044
diff changeset
851 octave_idx_type dist = 1;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
852 for (int i = 0; i < rank; i++)
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
853 dist *= dv(i);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
854
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
855 fftw_plan plan = octave_fftw_planner::create_plan (FFTW_FORWARD, rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
856 dv, 1, 1, dist, in, out);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
857
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
858 fftw_execute_dft (plan,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
859 reinterpret_cast<fftw_complex *> (const_cast<Complex *>(in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
860 reinterpret_cast<fftw_complex *> (out));
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
861
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
862 return 0;
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
863 }
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
864
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
865 int
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
866 octave_fftw::ifftNd (const Complex *in, Complex *out, const int rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
867 const dim_vector &dv)
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
868 {
5275
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe]
jwe
parents: 5044
diff changeset
869 octave_idx_type dist = 1;
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
870 for (int i = 0; i < rank; i++)
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
871 dist *= dv(i);
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
872
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
873 fftw_plan plan = octave_fftw_planner::create_plan (FFTW_BACKWARD, rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
874 dv, 1, 1, dist, in, out);
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
875
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
876 fftw_execute_dft (plan,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
877 reinterpret_cast<fftw_complex *> (const_cast<Complex *>(in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
878 reinterpret_cast<fftw_complex *> (out));
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
879
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
880 const size_t npts = dv.numel ();
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
881 const Complex scale = npts;
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
882 for (size_t i = 0; i < npts; i++)
4773
ccfbd6047a54 [project @ 2004-02-16 19:02:32 by jwe]
jwe
parents: 3874
diff changeset
883 out[i] /= scale;
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
884
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
885 return 0;
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
886 }
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
887
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
888 int
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
889 octave_fftw::fft (const float *in, FloatComplex *out, size_t npts,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
890 size_t nsamples, octave_idx_type stride, octave_idx_type dist)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
891 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
892 dist = (dist < 0 ? npts : dist);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
893
10350
12884915a8e4 merge MArray classes & improve Array interface
Jaroslav Hajek <highegg@gmail.com>
parents: 10314
diff changeset
894 dim_vector dv (npts, 1);
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
895 fftwf_plan plan = octave_float_fftw_planner::create_plan (1, dv, nsamples,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
896 stride, dist,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
897 in, out);
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
898
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
899 fftwf_execute_dft_r2c (plan, (const_cast<float *>(in)),
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
900 reinterpret_cast<fftwf_complex *> (out));
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
901
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
902 // Need to create other half of the transform.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
903
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
904 convert_packcomplex_1d (out, nsamples, npts, stride, dist);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
905
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
906 return 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
907 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
908
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
909 int
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
910 octave_fftw::fft (const FloatComplex *in, FloatComplex *out, size_t npts,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
911 size_t nsamples, octave_idx_type stride, octave_idx_type dist)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
912 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
913 dist = (dist < 0 ? npts : dist);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
914
10350
12884915a8e4 merge MArray classes & improve Array interface
Jaroslav Hajek <highegg@gmail.com>
parents: 10314
diff changeset
915 dim_vector dv (npts, 1);
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
916 fftwf_plan plan = octave_float_fftw_planner::create_plan (FFTW_FORWARD, 1,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
917 dv, nsamples,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
918 stride, dist,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
919 in, out);
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
920
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
921 fftwf_execute_dft (plan,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
922 reinterpret_cast<fftwf_complex *> (const_cast<FloatComplex *>(in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
923 reinterpret_cast<fftwf_complex *> (out));
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
924
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
925 return 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
926 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
927
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
928 int
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
929 octave_fftw::ifft (const FloatComplex *in, FloatComplex *out, size_t npts,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
930 size_t nsamples, octave_idx_type stride, octave_idx_type dist)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
931 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
932 dist = (dist < 0 ? npts : dist);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
933
10350
12884915a8e4 merge MArray classes & improve Array interface
Jaroslav Hajek <highegg@gmail.com>
parents: 10314
diff changeset
934 dim_vector dv (npts, 1);
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
935 fftwf_plan plan = octave_float_fftw_planner::create_plan (FFTW_BACKWARD, 1,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
936 dv, nsamples,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
937 stride, dist,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
938 in, out);
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
939
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
940 fftwf_execute_dft (plan,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
941 reinterpret_cast<fftwf_complex *> (const_cast<FloatComplex *>(in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
942 reinterpret_cast<fftwf_complex *> (out));
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
943
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
944 const FloatComplex scale = npts;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
945 for (size_t j = 0; j < nsamples; j++)
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
946 for (size_t i = 0; i < npts; i++)
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
947 out[i*stride + j*dist] /= scale;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
948
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
949 return 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
950 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
951
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
952 int
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
953 octave_fftw::fftNd (const float *in, FloatComplex *out, const int rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
954 const dim_vector &dv)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
955 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
956 octave_idx_type dist = 1;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
957 for (int i = 0; i < rank; i++)
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
958 dist *= dv(i);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
959
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
960 // Fool with the position of the start of the output matrix, so that
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
961 // creating other half of the matrix won't cause cache problems.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
962
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
963 octave_idx_type offset = (dv.numel () / dv(0)) * ((dv(0) - 1) / 2);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
964
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
965 fftwf_plan plan = octave_float_fftw_planner::create_plan (rank, dv, 1, 1,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
966 dist, in,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
967 out + offset);
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
968
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
969 fftwf_execute_dft_r2c (plan, (const_cast<float *>(in)),
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
970 reinterpret_cast<fftwf_complex *> (out+ offset));
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
971
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
972 // Need to create other half of the transform.
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
973
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
974 convert_packcomplex_Nd (out, dv);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
975
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
976 return 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
977 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
978
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
979 int
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
980 octave_fftw::fftNd (const FloatComplex *in, FloatComplex *out, const int rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
981 const dim_vector &dv)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
982 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
983 octave_idx_type dist = 1;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
984 for (int i = 0; i < rank; i++)
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
985 dist *= dv(i);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
986
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
987 fftwf_plan plan = octave_float_fftw_planner::create_plan (FFTW_FORWARD,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
988 rank, dv, 1, 1,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
989 dist, in, out);
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
990
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
991 fftwf_execute_dft (plan,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
992 reinterpret_cast<fftwf_complex *> (const_cast<FloatComplex *>(in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
993 reinterpret_cast<fftwf_complex *> (out));
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
994
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
995 return 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
996 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
997
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
998 int
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
999 octave_fftw::ifftNd (const FloatComplex *in, FloatComplex *out, const int rank,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
1000 const dim_vector &dv)
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1001 {
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1002 octave_idx_type dist = 1;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1003 for (int i = 0; i < rank; i++)
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1004 dist *= dv(i);
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1005
9516
fb933db0c517 convert fftw planner classes to singleton objects
John W. Eaton <jwe@octave.org>
parents: 8920
diff changeset
1006 fftwf_plan plan = octave_float_fftw_planner::create_plan (FFTW_BACKWARD,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
1007 rank, dv, 1, 1,
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
1008 dist, in, out);
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1009
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1010 fftwf_execute_dft (plan,
10314
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
1011 reinterpret_cast<fftwf_complex *> (const_cast<FloatComplex *>(in)),
07ebe522dac2 untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents: 10158
diff changeset
1012 reinterpret_cast<fftwf_complex *> (out));
7789
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1013
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1014 const size_t npts = dv.numel ();
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1015 const FloatComplex scale = npts;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1016 for (size_t i = 0; i < npts; i++)
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1017 out[i] /= scale;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1018
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1019 return 0;
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1020 }
82be108cc558 First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
1021
3828
adc217ebe692 [project @ 2001-05-02 14:50:46 by jwe]
jwe
parents:
diff changeset
1022 #endif