Skip to content

Commit

Permalink
gr-fft: Modernize fft code
Browse files Browse the repository at this point in the history 
* Add const where possible
* Disable copy assignment and copy constructor where not safe (for now) to copy
* Manual memory management -> smart pointers and vectors
* De-pointerify where possible
* assert -> static_assert
  • Loading branch information
@ThomasHabets @mbr0wn
ThomasHabets authored and mbr0wn committed Jan 5, 2020
1 parent ebe8ffc commit eae138c
Show file tree
Hide file tree
Showing 11 changed files with 90 additions and 122 deletions.
54 changes: 33 additions & 21 deletions gr-fft/include/gnuradio/fft/fft.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,13 +29,15 @@

#include <gnuradio/fft/api.h>
#include <gnuradio/gr_complex.h>
#include <volk/volk_alloc.hh>
#include <boost/thread.hpp>

namespace gr {
namespace fft {


/*! \brief Helper function for allocating complex* buffers
* TODO: Remove once the single user of this stops using it.
*/
FFT_API gr_complex* malloc_complex(int size);

Expand All @@ -48,6 +50,7 @@ FFT_API float* malloc_float(int size);
FFT_API double* malloc_double(int size);

/*! \brief Helper function for freeing fft buffers
* TODO: Remove once the single user of this stops using it.
*/
FFT_API void free(void* b);

Expand All @@ -71,26 +74,29 @@ class FFT_API planner
*/
class FFT_API fft_complex
{
int d_fft_size;
const int d_fft_size;
int d_nthreads;
gr_complex* d_inbuf;
gr_complex* d_outbuf;
volk::vector<gr_complex> d_inbuf;
volk::vector<gr_complex> d_outbuf;
void* d_plan;

public:
fft_complex(int fft_size, bool forward = true, int nthreads = 1);
// Copy disabled due to d_plan.
fft_complex(const fft_complex&) = delete;
fft_complex& operator=(const fft_complex&) = delete;
virtual ~fft_complex();

/*
* These return pointers to buffers owned by fft_impl_fft_complex
* into which input and output take place. It's done this way in
* order to ensure optimal alignment for SIMD instructions.
*/
gr_complex* get_inbuf() const { return d_inbuf; }
gr_complex* get_outbuf() const { return d_outbuf; }
gr_complex* get_inbuf() { return d_inbuf.data(); }
gr_complex* get_outbuf() { return d_outbuf.data(); }

int inbuf_length() const { return d_fft_size; }
int outbuf_length() const { return d_fft_size; }
int inbuf_length() const { return d_inbuf.size(); }
int outbuf_length() const { return d_outbuf.size(); }

/*!
* Set the number of threads to use for caclulation.
Expand All @@ -115,26 +121,29 @@ class FFT_API fft_complex
*/
class FFT_API fft_real_fwd
{
int d_fft_size;
const int d_fft_size;
int d_nthreads;
float* d_inbuf;
gr_complex* d_outbuf;
volk::vector<float> d_inbuf;
volk::vector<gr_complex> d_outbuf;
void* d_plan;

public:
fft_real_fwd(int fft_size, int nthreads = 1);
// Copy disabled due to d_plan.
fft_real_fwd(const fft_real_fwd&) = delete;
fft_real_fwd& operator=(const fft_real_fwd&) = delete;
virtual ~fft_real_fwd();

/*
* These return pointers to buffers owned by fft_impl_fft_real_fwd
* into which input and output take place. It's done this way in
* order to ensure optimal alignment for SIMD instructions.
*/
float* get_inbuf() const { return d_inbuf; }
gr_complex* get_outbuf() const { return d_outbuf; }
float* get_inbuf() { return d_inbuf.data(); }
gr_complex* get_outbuf() { return d_outbuf.data(); }

int inbuf_length() const { return d_fft_size; }
int outbuf_length() const { return d_fft_size / 2 + 1; }
int inbuf_length() const { return d_inbuf.size(); }
int outbuf_length() const { return d_outbuf.size(); }

/*!
* Set the number of threads to use for caclulation.
Expand All @@ -159,26 +168,29 @@ class FFT_API fft_real_fwd
*/
class FFT_API fft_real_rev
{
int d_fft_size;
const int d_fft_size;
int d_nthreads;
gr_complex* d_inbuf;
float* d_outbuf;
volk::vector<gr_complex> d_inbuf;
volk::vector<float> d_outbuf;
void* d_plan;

public:
fft_real_rev(int fft_size, int nthreads = 1);
// Copy disabled due to d_plan.
fft_real_rev(const fft_real_rev&) = delete;
fft_real_rev& operator=(const fft_real_rev&) = delete;
virtual ~fft_real_rev();

/*
* These return pointers to buffers owned by fft_impl_fft_real_rev
* into which input and output take place. It's done this way in
* order to ensure optimal alignment for SIMD instructions.
*/
gr_complex* get_inbuf() const { return d_inbuf; }
float* get_outbuf() const { return d_outbuf; }
gr_complex* get_inbuf() { return d_inbuf.data(); }
float* get_outbuf() { return d_outbuf.data(); }

int inbuf_length() const { return d_fft_size / 2 + 1; }
int outbuf_length() const { return d_fft_size; }
int inbuf_length() const { return d_inbuf.size(); }
int outbuf_length() const { return d_outbuf.size(); }

/*!
* Set the number of threads to use for caclulation.
Expand Down
2 changes: 1 addition & 1 deletion gr-fft/include/gnuradio/fft/goertzel.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,7 @@ class FFT_API goertzel
void set_params(int rate, int len, float freq);

// Process a input array
gr_complex batch(float* in);
gr_complex batch(const float* in);

// Process sample by sample
void input(const float& in);
Expand Down
4 changes: 2 additions & 2 deletions gr-fft/include/gnuradio/fft/goertzel_fc.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -45,9 +45,9 @@ class FFT_API goertzel_fc : virtual public sync_decimator

virtual void set_rate(int rate) = 0;

virtual float freq() = 0;
virtual float freq() const = 0;

virtual int rate() = 0;
virtual int rate() const = 0;
};

} /* namespace fft */
Expand Down
80 changes: 22 additions & 58 deletions gr-fft/lib/fft.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -47,7 +47,6 @@ static int my_fftw_read_char(void* f) { return fgetc((FILE*)f); }
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <cassert>
#include <stdexcept>

#include <boost/filesystem/operations.hpp>
Expand Down Expand Up @@ -167,37 +166,24 @@ static void export_wisdom()
// ----------------------------------------------------------------

fft_complex::fft_complex(int fft_size, bool forward, int nthreads)
: d_fft_size(fft_size), d_nthreads(nthreads), d_inbuf(fft_size), d_outbuf(fft_size)
{
// Hold global mutex during plan construction and destruction.
planner::scoped_lock lock(planner::mutex());

assert(sizeof(fftwf_complex) == sizeof(gr_complex));
static_assert(sizeof(fftwf_complex) == sizeof(gr_complex));

if (fft_size <= 0) {
throw std::out_of_range("fft_impl_fftw: invalid fft_size");
}

d_fft_size = fft_size;
d_inbuf = (gr_complex*)volk_malloc(sizeof(gr_complex) * inbuf_length(),
volk_get_alignment());
if (d_inbuf == 0) {
throw std::runtime_error("volk_malloc");
}
d_outbuf = (gr_complex*)volk_malloc(sizeof(gr_complex) * outbuf_length(),
volk_get_alignment());
if (d_outbuf == 0) {
volk_free(d_inbuf);
throw std::runtime_error("volk_malloc");
}

d_nthreads = nthreads;
config_threading(nthreads);
lock_wisdom();
import_wisdom(); // load prior wisdom from disk

d_plan = fftwf_plan_dft_1d(fft_size,
reinterpret_cast<fftwf_complex*>(d_inbuf),
reinterpret_cast<fftwf_complex*>(d_outbuf),
reinterpret_cast<fftwf_complex*>(d_inbuf.data()),
reinterpret_cast<fftwf_complex*>(d_outbuf.data()),
forward ? FFTW_FORWARD : FFTW_BACKWARD,
FFTW_MEASURE);

Expand All @@ -215,8 +201,6 @@ fft_complex::~fft_complex()
planner::scoped_lock lock(planner::mutex());

fftwf_destroy_plan((fftwf_plan)d_plan);
volk_free(d_inbuf);
volk_free(d_outbuf);
}

void fft_complex::set_nthreads(int n)
Expand All @@ -236,36 +220,28 @@ void fft_complex::execute() { fftwf_execute((fftwf_plan)d_plan); }
// ----------------------------------------------------------------

fft_real_fwd::fft_real_fwd(int fft_size, int nthreads)
: d_fft_size(fft_size),
d_nthreads(nthreads),
d_inbuf(fft_size),
d_outbuf(fft_size / 2 + 1)
{
// Hold global mutex during plan construction and destruction.
planner::scoped_lock lock(planner::mutex());

assert(sizeof(fftwf_complex) == sizeof(gr_complex));
static_assert(sizeof(fftwf_complex) == sizeof(gr_complex));

if (fft_size <= 0) {
throw std::out_of_range("gr::fft: invalid fft_size");
}

d_fft_size = fft_size;
d_inbuf = (float*)volk_malloc(sizeof(float) * inbuf_length(), volk_get_alignment());
if (d_inbuf == 0) {
throw std::runtime_error("volk_malloc");
}

d_outbuf = (gr_complex*)volk_malloc(sizeof(gr_complex) * outbuf_length(),
volk_get_alignment());
if (d_outbuf == 0) {
volk_free(d_inbuf);
throw std::runtime_error("volk_malloc");
}

d_nthreads = nthreads;
config_threading(nthreads);
lock_wisdom();
import_wisdom(); // load prior wisdom from disk

d_plan = fftwf_plan_dft_r2c_1d(
fft_size, d_inbuf, reinterpret_cast<fftwf_complex*>(d_outbuf), FFTW_MEASURE);
d_plan = fftwf_plan_dft_r2c_1d(fft_size,
d_inbuf.data(),
reinterpret_cast<fftwf_complex*>(d_outbuf.data()),
FFTW_MEASURE);

if (d_plan == NULL) {
fprintf(stderr, "gr::fft::fft_real_fwd: error creating plan\n");
Expand All @@ -281,8 +257,6 @@ fft_real_fwd::~fft_real_fwd()
planner::scoped_lock lock(planner::mutex());

fftwf_destroy_plan((fftwf_plan)d_plan);
volk_free(d_inbuf);
volk_free(d_outbuf);
}

void fft_real_fwd::set_nthreads(int n)
Expand All @@ -303,39 +277,31 @@ void fft_real_fwd::execute() { fftwf_execute((fftwf_plan)d_plan); }
// ----------------------------------------------------------------

fft_real_rev::fft_real_rev(int fft_size, int nthreads)
: d_fft_size(fft_size),
d_nthreads(nthreads),
d_inbuf(fft_size / 2 + 1),
d_outbuf(fft_size)
{
// Hold global mutex during plan construction and destruction.
planner::scoped_lock lock(planner::mutex());

assert(sizeof(fftwf_complex) == sizeof(gr_complex));
static_assert(sizeof(fftwf_complex) == sizeof(gr_complex));

if (fft_size <= 0) {
throw std::out_of_range("gr::fft::fft_real_rev: invalid fft_size");
}

d_fft_size = fft_size;
d_inbuf = (gr_complex*)volk_malloc(sizeof(gr_complex) * inbuf_length(),
volk_get_alignment());
if (d_inbuf == 0) {
throw std::runtime_error("volk_malloc");
}

d_outbuf = (float*)volk_malloc(sizeof(float) * outbuf_length(), volk_get_alignment());
if (d_outbuf == 0) {
volk_free(d_inbuf);
throw std::runtime_error("volk_malloc");
}

d_nthreads = nthreads;
config_threading(nthreads);
lock_wisdom();
import_wisdom(); // load prior wisdom from disk

// FIXME If there's ever a chance that the planning functions
// will be called in multiple threads, we've got to ensure single
// threaded access. They are not thread-safe.
d_plan = fftwf_plan_dft_c2r_1d(
fft_size, reinterpret_cast<fftwf_complex*>(d_inbuf), d_outbuf, FFTW_MEASURE);
d_plan = fftwf_plan_dft_c2r_1d(fft_size,
reinterpret_cast<fftwf_complex*>(d_inbuf.data()),
d_outbuf.data(),
FFTW_MEASURE);

if (d_plan == NULL) {
fprintf(stderr, "gr::fft::fft_real_rev: error creating plan\n");
Expand All @@ -351,8 +317,6 @@ fft_real_rev::~fft_real_rev()
planner::scoped_lock lock(planner::mutex());

fftwf_destroy_plan((fftwf_plan)d_plan);
volk_free(d_inbuf);
volk_free(d_outbuf);
}

void fft_real_rev::set_nthreads(int n)
Expand Down
22 changes: 10 additions & 12 deletions gr-fft/lib/fft_vcc_fftw.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -52,18 +52,16 @@ fft_vcc_fftw::fft_vcc_fftw(int fft_size,
io_signature::make(1, 1, fft_size * sizeof(gr_complex))),
d_fft_size(fft_size),
d_forward(forward),
d_fft(fft_size, forward, nthreads),
d_shift(shift)
{
d_fft = new fft_complex(d_fft_size, forward, nthreads);
if (!set_window(window))
throw std::runtime_error("fft_vcc: window not the same length as fft_size");
}

fft_vcc_fftw::~fft_vcc_fftw() { delete d_fft; }
void fft_vcc_fftw::set_nthreads(int n) { d_fft.set_nthreads(n); }

void fft_vcc_fftw::set_nthreads(int n) { d_fft->set_nthreads(n); }

int fft_vcc_fftw::nthreads() const { return d_fft->nthreads(); }
int fft_vcc_fftw::nthreads() const { return d_fft.nthreads(); }

bool fft_vcc_fftw::set_window(const std::vector<float>& window)
{
Expand All @@ -90,7 +88,7 @@ int fft_vcc_fftw::work(int noutput_items,

// copy input into optimally aligned buffer
if (!d_window.empty()) {
gr_complex* dst = d_fft->get_inbuf();
gr_complex* dst = d_fft.get_inbuf();
if (!d_forward && d_shift) {
unsigned int offset = (!d_forward && d_shift) ? (d_fft_size / 2) : 0;
int fft_m_offset = d_fft_size - offset;
Expand All @@ -103,30 +101,30 @@ int fft_vcc_fftw::work(int noutput_items,
}
} else {
if (!d_forward && d_shift) { // apply an ifft shift on the data
gr_complex* dst = d_fft->get_inbuf();
gr_complex* dst = d_fft.get_inbuf();
unsigned int len = (unsigned int)(floor(
d_fft_size / 2.0)); // half length of complex array
memcpy(&dst[0], &in[len], sizeof(gr_complex) * (d_fft_size - len));
memcpy(&dst[d_fft_size - len], &in[0], sizeof(gr_complex) * len);
} else {
memcpy(d_fft->get_inbuf(), in, input_data_size);
memcpy(d_fft.get_inbuf(), in, input_data_size);
}
}

// compute the fft
d_fft->execute();
d_fft.execute();

// copy result to our output
if (d_forward && d_shift) { // apply a fft shift on the data
unsigned int len = (unsigned int)(ceil(d_fft_size / 2.0));
memcpy(&out[0],
&d_fft->get_outbuf()[len],
&d_fft.get_outbuf()[len],
sizeof(gr_complex) * (d_fft_size - len));
memcpy(&out[d_fft_size - len],
&d_fft->get_outbuf()[0],
&d_fft.get_outbuf()[0],
sizeof(gr_complex) * len);
} else {
memcpy(out, d_fft->get_outbuf(), output_data_size);
memcpy(out, d_fft.get_outbuf(), output_data_size);
}

in += d_fft_size;
Expand Down

6 comments on commit eae138c

@yakkonaut
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@ThomasHabets I'm attempting to fix gr-iridium to work with gnuradio 3.9, and I want to make sure I understand what fft.h line 86 (fft_complex(const fft_complex&) = delete;) accomplishes? I think I can gut the calls to that out of the gr-iridium code but I want to make sure I communicate the reason for the change in the commit.

@ThomasHabets
Copy link
Contributor Author

@ThomasHabets ThomasHabets commented on eae138c Jul 1, 2020 •
edited

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@yakkonaut it means that fft_complex objects are not copy constructable. They're also not copy assignable (the line below).

This is because those objects have a member function d_plan, which is a raw pointer. There are then two options:

  1. Deep copy d_plan
  2. Disallow copies of fft_complex

Since d_plan is likely a somewhat expensive operation, and it's easy in C++ to accidentally create copies when you meant to move, or RVO, or something else, it's not only simpler but also safer to not allow copying.

If you copied fft_complex objects before this PR was merged, then that code was broken in that it destroyed d_plan twice if there was a copy, and gave a dangling pointer as soon as the first one was destroyed.

Looking at the code it seems you should just be able to remove the needless copy:

Change d_cfo_est_fft(fft::fft_complex(d_cfo_est_fft_size * d_fft_over_size_facor, true, 1)) to d_cfo_est_fft(d_cfo_est_fft_size * d_fft_over_size_facor, true, 1). I'm guessing that's the line that gives you the error?

@yakkonaut
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@ThomasHabets thank you so so so so much for looking at this. I didn't realize the fix was that simple... glad I stopped to ask for help. That is the line that's giving me the error... will keep pressing on with that change and make sure the tests still pass. Also thank you so much for explaining what's going on with d_plan there. It makes sense to me now why there had to be a change.

@yakkonaut
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There's one more issue in that file at line 278

// Temporary FFT to pre transform the filters
      fft::fft_complex fft_engine = fft::fft_complex(d_corr_fft_size);
      memset(fft_engine.get_inbuf(), 0, sizeof(gr_complex) * d_corr_fft_size);

It's not immediately apparent to me how to solve, but flagged this for the gr-iridium maintainers as well.

@ThomasHabets
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That one also looks needlessly verbose. Should just be fft::fft_complex fft_engine(d_corr_fft_size);.

@ThomasHabets
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

And you're welcome. :-)

Please sign in to comment.