Limiting to Vectors and Matrices necessary?

[jusack]

Hi @andreasnoack, thanks for your work on the rfft and pinging me in MagneticParticleImaging/MPIFiles.jl#136.

However, the limitation to 1D and 2D arrays in 0ae73f7 breaks our workflow and limits the applications of FFTA (from my perspective) unnecessarily. We use it to calculate the 1D FFT along a single dimension in a 4D array.

Would it be possible to revert this change? This would also allow the tests marked with @test_broken in the same commit to work again.

[andreasnoack]

Hi @jusack

The old implementation was not sufficiently tested and produced wrong results, see #76 (comment), so I decided to prioritize fixing the bug for a more narrow interface, see #76 (comment). These changes could be loosened, but

1. I think that the new restriction should only be relaxed when proper test cases are added.
2. I'd actually rather avoid this kind of implicit vectorization, see Autovectorization and AbstractFFTs #77 and specifically the examples in Autovectorization and AbstractFFTs #77 (comment), which might allow you to do what you want with the current version.

[jusack]

Ah, if its broken that makes sense to remove. And thank you for the example with mapslices. It feels a bit verbose, if you are used to the AbstractFFTs/FFTW way, but one can definitely work with that.

However, while trying to understand the issue with the code I played around a bit and noticed that apparently the region keyword for the rfft along one dimension of a 2D matrix also does not really work. The correct transform is calculated, but the result is the same as for the fft including the negative frequencies.

A = randn(5,5)
rfft(A) # 2D-transform, size (3,5) as expected
rfft(A,1) # 1D-transform along dimension 1, size (5,5) -> instead of (3,5)
rfft(A,2) # 1D-transform along dimension 2, size (5,5) -> instead of (5,3)
rfft(A,1)==fft(A,1) # true
rfft(A,1) == mapslices(rfft,A,dims=1) # false

This error was apparently also there before the recent changes. It seems like currently an assertion of N==NFFT seems to be necessary, until the underlying issue is fixed.

[andreasnoack]

I just realized this bug as well when adding error tests this afternoon. I believe I broke the 1D rfft on 2D arrays when fixing the 2D rfft. We need to decide on how to move forward. I'll start with getting all the error tests merged and then we should decide if we want to try to recover the FFTW "auto-vectorization" or delegate that to mapslices and restrict the code a bit more.

[jusack]

Oh, yeah you are right. It did not load the correct code when I though I switched to the previous version, it appears to be a new error. I would suggest trying to recover most of the auto-vectorization features for now, and think about removing them if it is decided to drop the AbstractFFTs dependency. From my side we can close this issue and continue the discussion in #77

[jusack]

I identified the mistake I made during the implementation in #69. I will implement a fix to try and restore the auto-vectorization for now.

[andreasnoack]

@jusack I'd actually started reinstantiating your version for the 1D rfft on ND arrays before seeing this, so if you haven't started working on this yet then you might want to hold off. I'll open a PR shortly.

[andreasnoack]

This was fixed by #86 and #87