@@ -14,39 +14,43 @@ RealFftPlan::RealFftPlan(int n)
1414}
1515
1616arr_cmplx RealFftPlan::solve (span_t <real_t > x) const {
17+ arr_cmplx r (n_);
18+ this ->solve (x, r);
19+ return r;
20+ }
21+
22+ void RealFftPlan::solve (span_t <real_t > x, mut_span_t <cmplx_t > r) const {
1723 using namespace std ::complex_literals;
1824 DSPLIB_ASSERT (x.size () == n_, " Input size must be equal FFT size" );
25+ // TODO: n/2+1 size support
26+ DSPLIB_ASSERT (r.size () == n_, " Output size must be equal FFT size" );
1927 const int n2 = n_ / 2 ;
2028
2129 arr_cmplx z (span (reinterpret_cast <const cmplx_t *>(x.data ()), n2));
2230 const auto Z = fft_->solve (z * 0.5 );
2331
24- arr_cmplx res (n_);
25-
2632 {
2733 const auto Xe = Z[0 ] + conj (Z[0 ]);
2834 const auto Xo = (conj (Z[0 ]) - Z[0 ]) * w_[0 ];
29- res [0 ].re = Xe.re - Xo.im ;
30- res [0 ].im = Xe.im + Xo.re ;
35+ r [0 ].re = Xe.re - Xo.im ;
36+ r [0 ].im = Xe.im + Xo.re ;
3137 }
3238
3339 for (int i = 1 ; i < n2; ++i) {
3440 const auto Zc = conj (Z[n2 - i]);
3541 const auto Xe = Z[i] + Zc;
3642 const auto Xo = (Zc - Z[i]) * w_[i];
37- res [i].re = Xe.re - Xo.im ;
38- res [i].im = Xe.im + Xo.re ;
39- res [n_ - i].re = res [i].re ;
40- res [n_ - i].im = -res [i].im ;
43+ r [i].re = Xe.re - Xo.im ;
44+ r [i].im = Xe.im + Xo.re ;
45+ r [n_ - i].re = r [i].re ;
46+ r [n_ - i].im = -r [i].im ;
4147 }
4248
4349 {
4450 const auto Xe = Z[0 ] + conj (Z[0 ]);
4551 const auto Xo = conj (Z[0 ]) - Z[0 ];
46- res [n2] = Xe.re + Xo.im ;
52+ r [n2] = Xe.re + Xo.im ;
4753 }
48-
49- return res;
5054}
5155
5256int RealFftPlan::size () const noexcept {
0 commit comments