xtensor-stack / xtensor-fftw Goto Github PK
View Code? Open in Web Editor NEWFFTW bindings for the xtensor C++14 multi-dimensional array library
License: BSD 3-Clause "New" or "Revised" License
FFTW bindings for the xtensor C++14 multi-dimensional array library
License: BSD 3-Clause "New" or "Revised" License
Greeting,
I wonder this package is not supporting the m1 chip now? Do you have any plans to implement?
Best,Fu
Hey @egpbos how about a new xtensor-fftw release?
I am happy to make the tag and update the recipes if you are ok with it.
Real FFTs on arrays where the last number in the shape is an odd number require special care. In particular: an inverse rFFT cannot automatically determine whether the output should be an odd or even sized (real) array, so the user has to supply this information. I added a flag for this in the relevant functions, but this behavior is undocumented. Add to docs.
Data arrays in FFTW must be in row major format (http://www.fftw.org/fftw3_doc/Row_002dmajor-Format.html). Currently, xtensor-fftw defines output arrays as row major, but it does not check for row major-ness on the input arrays. Should this be done or do we leave it to the user?
When following the install instructions: https://github.com/xtensor-stack/xtensor-fftw#installation
The required headers do not get installed, e.g. basic_float.hpp
, which is included from basic.hpp
does not get installed.
This is, I believe, because ${XTENSOR_FFTW_HEADERS}
is defined here without those headers:
https://github.com/xtensor-stack/xtensor-fftw/blob/master/CMakeLists.txt#L179
It's now included implicitly only if it happens to be in the same directory as xtensor... which it often does, but not necessarily. See xtensor CMakeLists.
Right now, due to my inexperience with templates, I have implemented the numpy-like interface functions using a number of (inline) functions that call a template instance, e.g.:
inline xt::xarray<std::complex<float> > fft (const xt::xarray<std::complex<float> > &input) {
return _fft_<std::complex<float>, std::complex<float>, 1, true, false, false, fftwf_plan_dft_1d, fftwf_execute, fftwf_destroy_plan> (input);
}
I'm not sure, but it seems likely that the compiler will simply compile all these functions when a user includes this header. This makes the resulting binary a lot larger than it needs to be and increases compile time.
I think I now know how to rewrite the _fft_
and _ifft_
templates in such a way that I can write using
template aliases for the three families (fft, rfft and hfft) that only depend on input precision. This will involve the following modifications to _fft_
and _ifft_
:
constexpr static
) in templates with parameters for precision, dimensionality, "c2c vs r2c" and direction.fftw_t
template as well (these are also constexpr static
).->
notation, using decltype
or declexpr
to determine the type based on template parameters for input precision, "c2c vs r2c" and direction.The latter point was the crucial knowledge gap that caused me to write the functions as I did, since I didn't get auto-type-deduction working for the output type.
Currently, the helper functions in helper.hpp
only work on single dimensional arrays. We should extend this to multi-dimensional arrays. numpy.fft
's helper functions use the axes
keyword argument for this. This syntax could also be used in broadcasting operations (#4). Possibly, some other syntax is more appropriate in the C++ / xtensor context. We should investigate.
A cling Jupyter notebook showing interactive usage. Something with an image, something with filtering, maybe something with audio. Maybe there are numpy.fft notebooks that we could draw inspiration from as well.
On my machine (Rocky Linux 8.9; FFTW 3.3.5), when loading this data and execute FFTW, it reports a smallbin double linked list corruption.
The example is as following (C++ code):
fstream fs_debug;
fs_debug.open("/home/winfred/l1b-pro/data/output/debug-8-7298.csv", ios::in);
xtensor<float, 2> delta_sigx_han_t = xt::load_csv<float>(fs_debug);
xtensor<float, 1> delta_sigx_han_tt = xt::zeros<float>({delta_sigx_han_t.size()});
for (int i = 0; i < delta_sigx_han_t.size(); i++)
{
delta_sigx_han_tt(i) = delta_sigx_han_t(0, i);
}
auto pdsig_t = xt::fftw::rfft(delta_sigx_han_tt);
FFTW is not thread safe, except for fftw_execute
.
How about using c++11 locks around calls to fftw_plan_dft
and fftw_destroy_plan
for thread safety?
Thank you for developing xtensor-fftw
.
I would like to request the following extension that will be improve its usability.
It will be better if there is a broadcasting in xtensor-fftw
.
I'm studying on multi-channel acoustic signal processing, and often have to conduct FFTs to microphones.
In numpy, we can broadcast FFT like this way
np.fft.rfft(np.ndarray([M, T]), axis=1).shape
In my case M and T indicates the number of microphones and time frames.
This will be simplify codes of such signal processings.
Best,
Take inspiration from the other xtensor pages, e.g. xensor-io has a nice compact one, that would be good as a header. The rest can be something like:
CMake fails to find FFTW (version 3.3.9) which was installed from source. I use CMake version 3.19.3 on Ubuntu 20.04.
$ cmake --version
cmake version 3.19.3
CMake suite maintained and supported by Kitware (kitware.com/cmake).
FFTW is installed in /usr/local
:
/usr/local/lib/libfftw3q.a
/usr/local/lib/libfftw3q.a
/usr/local/include/fftw3.h
/usr/local/include/fftw3.f
/usr/local/include/fftw3l.f03
/usr/local/include/fftw3q.f03
/usr/local/include/fftw3.f03
/usr/local/lib/pkgconfig/fftw3q.pc
/usr/local/lib/cmake/fftw3q/FFTW3qConfig.cmake
/usr/local/lib/cmake/fftw3q/FFTW3qConfigVersion.cmake
/usr/local/lib/cmake/fftw3q/FFTW3LibraryDepends.cmake
/usr/local/lib/cmake/fftw3q/FFTW3LibraryDepends-release.cmake
I suspect that your CMake code is no longer compatible with current CMake versions.
Simple cmake call:
$ cmake ../xtensor-fftw-0.2.6
-- The C compiler identification is GNU 9.3.0
-- The CXX compiler identification is GNU 9.3.0
-- Detecting C compiler ABI info
-- Detecting C compiler ABI info - done
-- Check for working C compiler: /usr/bin/cc - skipped
-- Detecting C compile features
-- Detecting C compile features - done
-- Detecting CXX compiler ABI info
-- Detecting CXX compiler ABI info - done
-- Check for working CXX compiler: /usr/bin/c++ - skipped
-- Detecting CXX compile features
-- Detecting CXX compile features - done
-- Building xtensor-fftw v0.2.6
-- Performing Test arch_native_supported
-- Performing Test arch_native_supported - Success
-- Found xtensor: /usr/local/include/xtensor
-- Found xtl: /usr/local/include/xtl
-- Found PkgConfig: /usr/bin/pkg-config (found version "0.29.1")
CMake Error at /snap/cmake/769/share/cmake-3.19/Modules/FindPackageHandleStandardArgs.cmake:218 (message):
Could NOT find FFTW (missing: FLOAT_LIB DOUBLE_LIB LONGDOUBLE_LIB)
Call Stack (most recent call first):
/snap/cmake/769/share/cmake-3.19/Modules/FindPackageHandleStandardArgs.cmake:582 (_FPHSA_FAILURE_MESSAGE)
cmake/Modules/findFFTW/FindFFTW.cmake:277 (find_package_handle_standard_args)
CMakeLists.txt:88 (find_package)
-- Configuring incomplete, errors occurred!
See also "/home/bebuch/media/tmp/install/xtensor-fftw-0.2.6-build/CMakeFiles/CMakeOutput.log".
Call with FFTW_ROOT:
$ cmake -DFFTW_ROOT=/usr/local ../xtensor-fftw-${XTENSOR_FFTW_VERSION}
-- The C compiler identification is GNU 9.3.0
-- The CXX compiler identification is GNU 9.3.0
-- Detecting C compiler ABI info
-- Detecting C compiler ABI info - done
-- Check for working C compiler: /usr/bin/cc - skipped
-- Detecting C compile features
-- Detecting C compile features - done
-- Detecting CXX compiler ABI info
-- Detecting CXX compiler ABI info - done
-- Check for working CXX compiler: /usr/bin/c++ - skipped
-- Detecting CXX compile features
-- Detecting CXX compile features - done
-- Building xtensor-fftw v0.2.6
-- Performing Test arch_native_supported
-- Performing Test arch_native_supported - Success
-- Found xtensor: /usr/local/include/xtensor
-- Found xtl: /usr/local/include/xtl
CMake Warning (dev) at CMakeLists.txt:88 (find_package):
Policy CMP0074 is not set: find_package uses <PackageName>_ROOT variables.
Run "cmake --help-policy CMP0074" for policy details. Use the cmake_policy
command to set the policy and suppress this warning.
CMake variable FFTW_ROOT is set to:
/usr/local
For compatibility, CMake is ignoring the variable.
This warning is for project developers. Use -Wno-dev to suppress it.
-- Found PkgConfig: /usr/bin/pkg-config (found version "0.29.1")
CMake Error at /snap/cmake/769/share/cmake-3.19/Modules/FindPackageHandleStandardArgs.cmake:218 (message):
Could NOT find FFTW (missing: FLOAT_LIB DOUBLE_LIB LONGDOUBLE_LIB)
Call Stack (most recent call first):
/snap/cmake/769/share/cmake-3.19/Modules/FindPackageHandleStandardArgs.cmake:582 (_FPHSA_FAILURE_MESSAGE)
cmake/Modules/findFFTW/FindFFTW.cmake:277 (find_package_handle_standard_args)
CMakeLists.txt:88 (find_package)
-- Configuring incomplete, errors occurred!
See also "/home/bebuch/media/tmp/install/xtensor-fftw-0.2.6-build/CMakeFiles/CMakeOutput.log".
We want to implement a very simple FFT interface, mimicking the numpy.fft default behavior. This means that we want the following matrix of functions:
fft
for standard complex to complex transforms, rfft
for real FFTs (real signal) and hfft
for Hermitian FFTs (real Fourier spectrum);i
, e.g. ifft
;2
, 3
or n
, e.g. irfft3
.That's a total of 24 functions.
Additionally, we want to support single, double and long double precision arrays. This can be implemented through template specializations. The total amount of implemented functions would then be 72.
Also part of this interface would be broadcasting over multiple array dimensions, if the dimension of the function is smaller than that of the array itself; see #4.
The test on my last PR failed in a way that does not match the PR. Added to that, the error is the same than the last failed run on master.
My PR failure:
Last build on master:
I may fix it if you dont before tomorrow (note Im in the US west coast so my schedule does not overlap with the main developer here)
We should benchmark whether the 1-3 dimensional special functions are actually faster than the n-dimensional one. If not, we can just use one implementation function for all ranks, with simple wrappers for the 1, 2 and 3 dimensional cases. This would reduce the number of implemented functions from 72 to 18, plus 54 simple wrapper functions.
Following #28, to make things even more general, replace shape type with xarray<T, L>::shape_type
when you get the shape of an array, thus the code is more flexible and less sensitive to changes of shape / strides types in xtensor (suggestion by @JohanMabille ).
Follow some good suggestions from https://youtu.be/bsXLMQ6WgIk.
add_library
for installation of header only libraries", which xtensor-fftw is (xtensor itself as well).CXX_CLANG_TIDY=ON
by default.Nevermind, misinterpreted some stuff.
As in https://github.com/QuantStack/xtensor-io/blob/a57e3b90ca1025a180c5813c56b8a344ec2d8165/include/xtensor-io/xaudio.hpp#L20-L21, we can put the library loading step in the headers themselves, which saves cling/notebook users from doing it.
After implementing the basic numpy.fft-like interface (#6), we should benchmark it against the actual numpy.fft functions.
We use template specializations to make sure we get the correct precision. The problem, for instance, with a non-template f(const xt::xarray<float> &input)
is that it will also compile when you pass a xt::xarray<double>
. Passing by const-reference in this sense behaves similarly to passing by value; it triggers the creation of a temporary variable -- input in this case -- though in the case of a reference the data is not actually copied.
Using = delete
syntax on the default template, e.g.:
template<typename real_t>
xt::xarray< std::complex<real_t> > fft(const xt::xarray<real_t> &input) = delete;
we can make sure that calls to non-implemented specializations don't compile. If this is left out, the compilation will succeed, but the linker will fail, and this gives less informative error messages.
Unfortunately, clang has a bug that prevents us from using the nice = delete
syntax. For the moment we're using a workaround:
template<typename real_t>
xt::xarray< std::complex<real_t> > fft(const xt::xarray<real_t> &input) {
static_assert(sizeof(real_t) == 0, "Only specializations of fft can be used");
};
This is ugly. Once the bug in clang is fixed, we will want to replace this with the = delete
syntax.
Look into which license xtensor-fftw is exactly supposed/allowed to use for what (code, conda package), given that FFTW3 itself is under GPL.
The simple numpy-style interface (#6) will by default hide some of FFTW's features, opting for simplicity and default settings for decent performance. We will also implement a high performance interface to access functions like
In basic.hpp, there are a few exceptions to safeguard against nullptr
s from FFTW plan creation, which for instance can happen if you try to build a plan for a multi-dimensional real FFT with a FFTW_PRESERVE_INPUT
flag, because such algorithms are not available in FFTW. These exceptions are there more as a developer safeguard than a user safeguard. Eventually, it would be nicer if such "illegal" combinations of plan parameters were disallowed by the code itself, e.g. by not allowing users to manually specify planner flags, but to hard code all legal possibilities. This will eliminate the need for the exceptions.
Therefore, note to self: don't try to get 100% coverage by making a test with illegal input, or something like that. Just fix the API so that the exceptions can be removed.
The current setup is to use Travis for Linux and macOS and Appveyor for Windows. We should replace both by just GitHub Actions to keep things simple. Also, Travis doesn't actually work anymore, so Actions would fix that as well.
Hi, I'm adding xtensor-fftw to our build process at work. We don't use FFTW, but use a compatible implementation through MKL. The CMakeLists file however doesn't check for this. Could you please add a check for MKL? Thanks!
Use e.g. licecap. Put it on Twitter and maybe on the documentation page under Interactive examples.
I am currently using xtensor
as well as xtensor-blas
. To use xtensor-fftw
, Conda
asks me to downgrade these packages:
Package Version Build Channel Size
───────────────────────────────────────────────────────────────
Install:
───────────────────────────────────────────────────────────────
+ xtensor-fftw 0.2.6 h4bd325d_1 conda-forge 12kB
Downgrade:
───────────────────────────────────────────────────────────────
- xtl 0.7.5 hf52228f_0 conda-forge Cached
+ xtl 0.6.23 h4bd325d_1 conda-forge 98kB
- xtensor 0.24.7 h00ab1b0_0 conda-forge Cached
+ xtensor 0.21.10 h4bd325d_0 conda-forge 182kB
- xtensor-blas 0.20.0 h4bd325d_0 conda-forge Cached
+ xtensor-blas 0.17.2 h998ca80_1 conda-forge 264kB
Summary:
Install: 1 packages
Downgrade: 3 packages
Is it possible to use xtensor-fftw
with recent xtensor
versions or are there some issues?
We should add items to the numpy to xtensor cheat sheet and a code snippet in the xtensor related projects section.
Some decent machine + maybe a website to publish the recorded results. Continuous benchmarks would be even better, since then we can catch performance regressions more easily.
Hi,
First of all, thank you for developing both xtensor
and xtensor-fftw
! I encounter an issue while trying to integrate xtensor-fftw
with my current project using xtensor
(for which I install both xtensor
and xtensor-fftw
through conda install -c conda-forge xtensor xsimd xtensor-fftw
) --
Here is my CMakeLists.txt
file --
cmake_minimum_required(VERSION 3.1)
project(example)
find_package(xtl REQUIRED)
find_package(xtensor REQUIRED)
find_package(xtensor-fftw REQUIRED)
set(CMAKE_CXX_STANDARD 17)
set(CMAKE_CXX_STANDARD_REQUIRED TRUE)
set(CMAKE_CXX_EXTENSIONS TRUE)
add_executable(main src/main.cpp)
target_link_libraries(main xtensor xtensor::optimize xtensor::use_xsimd)
And, the main.cpp
looks like this:
#include <complex>
#include <iostream>
#include <xtensor-fftw/basic.hpp>
#include <xtensor-fftw/helper.hpp>
#include <xtensor/xarray.hpp>
#include <xtensor/xbuilder.hpp>
#include <xtensor/xio.hpp>
#include <xtensor/xmath.hpp>
int main(int argc, char *argv[]) {
using T = std::complex<double>;
auto size = 128;
xt::xarray<T> arr = xt::zeros<T>({size, 1});
for (auto i = 0; i < 10; i++) {
arr(i, 0) = arr(size - i, 0) = 1.0;
}
auto arrfft = xt::fftw::fft(arr);
std::cout << arrfft;
return 0;
}
However, the linker
complaints
✖ make (xt)
Consolidate compiler generated dependencies of target main
[ 50%] Building CXX object CMakeFiles/main.dir/src/main.cpp.o
[100%] Linking CXX executable main
Undefined symbols for architecture x86_64:
"_fftw_destroy_plan", referenced from:
xt::xarray_container<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >, (xt::layout_type)1, xt::svector<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type, 4ul, std::__1::allocator<std::__1::vector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type>, true>, xt::xtensor_expression_tag> xt::fftw::_fft_<std::__1::complex<double>, std::__1::complex<double>, 1ul, -1, true, false, false, fftw_plan_dft_1d, fftw_execute, fftw_destroy_plan, void>(xt::xarray_container<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >, (xt::layout_type)1, xt::svector<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type, 4ul, std::__1::allocator<std::__1::vector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type>, true>, xt::xtensor_expression_tag> const&) in main.cpp.o
"_fftw_execute", referenced from:
xt::xarray_container<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >, (xt::layout_type)1, xt::svector<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type, 4ul, std::__1::allocator<std::__1::vector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type>, true>, xt::xtensor_expression_tag> xt::fftw::_fft_<std::__1::complex<double>, std::__1::complex<double>, 1ul, -1, true, false, false, fftw_plan_dft_1d, fftw_execute, fftw_destroy_plan, void>(xt::xarray_container<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >, (xt::layout_type)1, xt::svector<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type, 4ul, std::__1::allocator<std::__1::vector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type>, true>, xt::xtensor_expression_tag> const&) in main.cpp.o
"_fftw_plan_dft_1d", referenced from:
std::__1::enable_if<(dimensional::is_1<1ul, true>::value) && ((-1) != (0)), xt::fftw::fftw_t<std::__1::complex<double> >::plan>::type xt::fftw::fftw_plan_dft_caller<1ul, -1, true, std::__1::complex<double>, std::__1::complex<double>, fftw_plan_dft_1d, false, false>(xt::xarray_container<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >, (xt::layout_type)1, xt::svector<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type, 4ul, std::__1::allocator<std::__1::vector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type>, true>, xt::xtensor_expression_tag> const&, xt::xarray_container<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >, (xt::layout_type)1, xt::svector<xt::uvector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type, 4ul, std::__1::allocator<std::__1::vector<std::__1::complex<double>, xsimd::aligned_allocator<std::__1::complex<double>, 32ul> >::size_type>, true>, xt::xtensor_expression_tag>&, unsigned int, bool) in main.cpp.o
ld: symbol(s) not found for architecture x86_64
clang: error: linker command failed with exit code 1 (use -v to see invocation)
make[2]: *** [main] Error 1
make[1]: *** [CMakeFiles/main.dir/all] Error 2
make: *** [all] Error 2
I'm wondering could you please help me with this?
The CMakeFiles have become a bit messy and seem overcomplicated. I think a modern CMake makeover would help a lot. Things like FFTW_LINK_FLAGS
, for instance, should not be used as such, since linker flags and linking targets are two things that are defined separately, so combining them in one variable is not CMake-ish. Other things like set
for custom variables are also code smells.
The default behavior for inverse real FFT functions (c2r) in FFTW is to destroy the input arrays during the calculation. This is probably very unexpected for most users.
We should add this to the documentation, probably multiple times (in the general howtos, and in the function docs). Maybe we should even print run-time warnings.
I don't know if it makes sense and if it's easy to achieve it, but it would be nice to not have to evaluate the xarrays for a call of fft function, meaning that the implementation of the functions would be:
template <class T>
inline xt::xarray<std::complex<float> > fft3 (const xt::xexpression<T> &input) {
return ...;
}
instead of
inline xt::xarray<std::complex<float> > fft3 (const xt::xarray<std::complex<float> > &input) {
return ...;
}
And maybe the return type could be an xexpression too..
A declarative, efficient, and flexible JavaScript library for building user interfaces.
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google ❤️ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.