Skip to content

Commit 62dfa6f

Browse files
vitalsongvitalsong
committed
FFT span optimization
1 parent 9034733 commit 62dfa6f

File tree

7 files changed

+76
-52
lines changed

7 files changed

+76
-52
lines changed

‎lib/fft/cmplx-ifft.h‎

Lines changed: 19 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -12,26 +12,32 @@ class CmplxIfftPlan : public IfftPlanC
1212
}
1313

1414
arr_cmplx solve(span_t<cmplx_t> x) const final {
15-
const real_t m = real_t(1) / x.size();
16-
arr_cmplx y(x);
17-
y *= m;
18-
_inplace_conj(y);
19-
y = fft_->solve(y);
20-
_inplace_conj(y);
21-
return y;
15+
arr_cmplx r(x);
16+
this->solve(x, r);
17+
return r;
18+
}
19+
20+
void solve(span_t<cmplx_t> x, mut_span_t<cmplx_t> r) const final {
21+
const int n = fft_->size();
22+
DSPLIB_ASSERT(x.size() == n, "array size error");
23+
DSPLIB_ASSERT(x.size() == r.size(), "array size error");
24+
arr_cmplx t(n);
25+
const real_t m = real_t(1) / n;
26+
for (int i = 0; i < n; ++i) {
27+
t[i].re = x[i].re * m;
28+
t[i].im = -(x[i].im * m);
29+
}
30+
fft_->solve(t, r);
31+
for (int i = 0; i < n; ++i) {
32+
r[i].im = -r[i].im;
33+
}
2234
}
2335

2436
int size() const noexcept final {
2537
return fft_->size();
2638
}
2739

2840
private:
29-
static void _inplace_conj(arr_cmplx& x) {
30-
for (auto& v : x) {
31-
v.im = -v.im;
32-
}
33-
}
34-
3541
std::shared_ptr<FftPlanC> fft_;
3642
};
3743

‎lib/fft/fact-fft.cpp‎

Lines changed: 16 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -59,12 +59,12 @@ class PlanTree
5959
return _n;
6060
}
6161

62-
[[nodiscard]] PlanTree* q_plan() const noexcept {
62+
[[nodiscard]] const PlanTree* q_plan() const noexcept {
6363
assert(has_next());
6464
return _q;
6565
}
6666

67-
[[nodiscard]] PlanTree* p_plan() const noexcept {
67+
[[nodiscard]] const PlanTree* p_plan() const noexcept {
6868
assert(has_next());
6969
return _p;
7070
}
@@ -102,8 +102,8 @@ class PlanTree
102102
}
103103

104104
const int _n;
105-
PlanTree* _p{nullptr};
106-
PlanTree* _q{nullptr};
105+
const PlanTree* _p{nullptr};
106+
const PlanTree* _q{nullptr};
107107
std::shared_ptr<FftPlanC> _solver;
108108
};
109109

@@ -175,19 +175,26 @@ void _facfft(const PlanTree* plan, cmplx_t* restrict x, cmplx_t* restrict mem, c
175175
//-----------------------------------------------------------------------------------------------------------------------------
176176
FactorFFTPlan::FactorFFTPlan(int n)
177177
: _n{n}
178-
, _px(n) {
178+
, _twiddle{expj(-2 * pi * arange(n) / n)} //TODO: only part of the table is needed
179+
{
179180
DSPLIB_ASSERT(!isprime(n), "fft size must not be a prime number");
180-
_twiddle = expj(-2 * pi * arange(n) / n); //TODO: only part of the table is needed
181181
_plan = std::make_shared<PlanTree>(n);
182182
}
183183

184184
[[nodiscard]] arr_cmplx FactorFFTPlan::solve(span_t<cmplx_t> x) const {
185-
DSPLIB_ASSERT(x.size() == _n, "input vector size is not equal fft size");
186-
arr_cmplx r(x); //TODO: remove copy
187-
_facfft(_plan.get(), r.data(), _px.data(), _twiddle.data(), _n);
185+
arr_cmplx r(_n);
186+
this->solve(x, r);
188187
return r;
189188
}
190189

190+
void FactorFFTPlan::solve(span_t<cmplx_t> x, mut_span_t<cmplx_t> r) const {
191+
DSPLIB_ASSERT(x.size() == _n, "input array size is not equal fft size");
192+
DSPLIB_ASSERT(x.size() == r.size(), "output array size error");
193+
arr_cmplx tmp(_n);
194+
r = x; //TODO: remove copy
195+
_facfft(_plan.get(), r.data(), tmp.data(), _twiddle.data(), _n);
196+
}
197+
191198
[[nodiscard]] int FactorFFTPlan::size() const noexcept {
192199
return _n;
193200
}

‎lib/fft/fact-fft.h‎

Lines changed: 4 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -2,7 +2,6 @@
22

33
#include <dsplib/math.h>
44
#include <dsplib/fft.h>
5-
#include <dsplib/math.h>
65

76
#include <memory>
87

@@ -15,14 +14,14 @@ class FactorFFTPlan : public FftPlanC
1514
{
1615
public:
1716
explicit FactorFFTPlan(int n);
18-
~FactorFFTPlan() = default;
17+
~FactorFFTPlan() override = default;
1918
[[nodiscard]] arr_cmplx solve(span_t<cmplx_t> x) const final;
19+
void solve(span_t<cmplx_t> x, mut_span_t<cmplx_t> r) const final;
2020
[[nodiscard]] int size() const noexcept final;
2121

2222
private:
23-
int _n;
24-
arr_cmplx _twiddle;
25-
mutable arr_cmplx _px; ///< tmp matrix for transpose
23+
const int _n;
24+
const arr_cmplx _twiddle;
2625
std::shared_ptr<PlanTree> _plan;
2726
};
2827

‎lib/fft/pow2-fft.cpp‎

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -65,10 +65,10 @@ void _bitreverse(const cmplx_t* restrict x, cmplx_t* restrict y, const int32_t*
6565

6666
Pow2FftPlan::Pow2FftPlan(int n)
6767
: n_{n}
68-
, l_{nextpow2(n_)} {
68+
, l_{nextpow2(n_)}
69+
, bitrev_{_gen_bitrev_table(n)}
70+
, coeffs_{_gen_coeffs_table(n)} {
6971
DSPLIB_ASSERT(ispow2(n), "FFT size must be power of 2");
70-
bitrev_ = _gen_bitrev_table(n);
71-
coeffs_ = _gen_coeffs_table(n);
7272
}
7373

7474
void Pow2FftPlan::solve(span_t<cmplx_t> x, mut_span_t<cmplx_t> r) const {

‎lib/fft/pow2-fft.h‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -19,8 +19,8 @@ class Pow2FftPlan : public FftPlanC
1919

2020
const int n_;
2121
const int l_;
22-
std::vector<int32_t> bitrev_;
23-
std::vector<cmplx_t> coeffs_;
22+
const std::vector<int32_t> bitrev_;
23+
const std::vector<cmplx_t> coeffs_;
2424
};
2525

2626
} // namespace dsplib

‎lib/fft/small-fft.h‎

Lines changed: 27 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -1,7 +1,5 @@
11
#pragma once
22

3-
#include "dsplib/array.h"
4-
#include "dsplib/types.h"
53
#include <dsplib/fft.h>
64
#include <dsplib/assert.h>
75

@@ -22,30 +20,35 @@ class SmallFftC : public FftPlanC
2220
~SmallFftC() override {
2321
}
2422

25-
[[nodiscard]] arr_cmplx solve(span_t<cmplx_t> x) const final {
23+
void solve(span_t<cmplx_t> x, mut_span_t<cmplx_t> r) const final {
2624
DSPLIB_ASSERT(x.size() == n_, "input size error");
27-
arr_cmplx y(x.size());
25+
DSPLIB_ASSERT(x.size() == r.size(), "input size error");
2826
switch (n_) {
2927
case 1:
30-
y[0] = x[0];
28+
r[0] = x[0];
3129
break;
3230
case 2:
33-
_fft_n2(x.data(), y.data());
31+
_fft_n2(x.data(), r.data());
3432
break;
3533
case 3:
36-
_fft_n3(x.data(), y.data());
34+
_fft_n3(x.data(), r.data());
3735
break;
3836
case 4:
39-
_fft_n4(x.data(), y.data());
37+
_fft_n4(x.data(), r.data());
4038
break;
4139
case 8:
42-
_fft_n8(x.data(), y.data());
40+
_fft_n8(x.data(), r.data());
4341
break;
4442
default:
4543
DSPLIB_THROW("size not supported");
4644
break;
4745
}
48-
return y;
46+
}
47+
48+
[[nodiscard]] arr_cmplx solve(span_t<cmplx_t> x) const final {
49+
arr_cmplx r(x.size());
50+
this->solve(x, r);
51+
return r;
4952
}
5053

5154
[[nodiscard]] int size() const noexcept final {
@@ -137,30 +140,36 @@ class SmallFftR : public FftPlanR
137140
~SmallFftR() override {
138141
}
139142

140-
[[nodiscard]] arr_cmplx solve(span_t<real_t> x) const final {
143+
void solve(span_t<real_t> x, mut_span_t<cmplx_t> r) const final {
141144
DSPLIB_ASSERT(x.size() == n_, "input size error");
142-
arr_cmplx y(x.size());
145+
DSPLIB_ASSERT(x.size() == r.size(), "input size error");
143146
switch (n_) {
144147
case 1:
145-
y[0] = x[0];
148+
r[0].re = x[0];
149+
r[0].im = 0;
146150
break;
147151
case 2:
148-
_fft_n2(x.data(), y.data());
152+
_fft_n2(x.data(), r.data());
149153
break;
150154
case 3:
151-
_fft_n3(x.data(), y.data());
155+
_fft_n3(x.data(), r.data());
152156
break;
153157
case 4:
154-
_fft_n4(x.data(), y.data());
158+
_fft_n4(x.data(), r.data());
155159
break;
156160
case 8:
157-
_fft_n8(x.data(), y.data());
161+
_fft_n8(x.data(), r.data());
158162
break;
159163
default:
160164
DSPLIB_THROW("size not supported");
161165
break;
162166
}
163-
return y;
167+
}
168+
169+
[[nodiscard]] arr_cmplx solve(span_t<real_t> x) const final {
170+
arr_cmplx r(x.size());
171+
this->solve(x, r);
172+
return r;
164173
}
165174

166175
[[nodiscard]] int size() const noexcept final {

‎lib/subband.cpp‎

Lines changed: 5 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -182,11 +182,12 @@ class ChannelizerImpl : public DFTFilterBank
182182
const auto* restrict buf = buf_[decim_ * k];
183183
const auto* restrict flt = fview_[k];
184184
for (int m = 0; m < nbands_; m++) {
185-
pout[nbands_ - m - 1] += flt[m] * buf[m];
185+
pout[m] += flt[m] * buf[m];
186186
}
187187
}
188188

189-
return fft_->solve(pout);
189+
//TODO: remove conj
190+
return conj(fft_->solve(pout));
190191
}
191192

192193
private:
@@ -212,6 +213,8 @@ class ChannelSynthesizerImpl : public DFTFilterBank
212213

213214
auto xx = ifft_->solve(x);
214215
xx *= nbands_;
216+
217+
//TODO: fft and flip?
215218
buf_.push(xx, true);
216219

217220
// calculate outputs of polyphase filters

0 commit comments

Comments
 (0)