[Octave-bug-tracker] [bug #49026] function fft return wrong result or crash if n=1 and dim>2

[Rik]

Update of bug #49026 (project octave):

                  Status:               Confirmed => Patch Submitted        

    _______________________________________________________

Follow-up Comment #2:

The problem seems to be in do_fft in libinterp/corefcn/fft.cc.  I quote the
code, along for the first case for real FFT of single matrices.


  if (n_points < 0)
    n_points = dims(dim);
  else
    dims(dim) = n_points;

  if (dims.any_zero () || n_points == 0)
    {
      if (arg.is_single_type ())
        return octave_value (FloatNDArray (dims));
      else
        return octave_value (NDArray (dims));
    }

  if (arg.is_single_type ())
    {
      if (arg.is_real_type ())
        {
          FloatNDArray nda = arg.float_array_value ();

          nda.resize (dims, 0.0);
          retval = (type != 0 ? nda.ifourier (dim) : nda.fourier (dim));
        }


The original array was a 4-D array (2x4x2x2).  The call to fft


fft (x,1,4)


has asked to use just a single point along the fourth (last) dimension.  In
the C++ code,


    dims(dim) = n_points;


the final dimension is set to 1 so the new dimensions are 2x4x2x1.

Eventually, nda.resize is called, and Octave knows to trim singleton
dimensions, so the matrix becomes a 3-D matrix 2x4x2.  You can check and find
that fft works correctly when the dimension is not the trailing dimension.

That's the problem.  I can see a few solutions.  First, from someone who
understands the FFT well, if I take an FFT with a single point, is the result
really just an indexing operation?  In this example


fft (x,1,4) === x(:,:,:,1)  ???


If that is the case, then maybe we can implement that equivalent index
operation on an octave_value in C++ only when n_points == 1.  This wouldn't be
too much of a hassle because there is already a special case for n_points ==
0.

Second, maybe there is a way to resize an array and avoid trimming singleton
dimensions.

Third, and least interesting, each fourier and ifourier routine in liboctave
needs to be updated to handle this special case.

I'm attaching a patch that works for me.  Mike, can you test it? 

Since this is a crash I would apply it to the stable branch if it passes QC.

(file #38789)
    _______________________________________________________

Additional Item Attachment:

File name: fft.cset                       Size:4 KB


    _______________________________________________________

Reply to this item at:

  <http://savannah.gnu.org/bugs/?49026>

_______________________________________________
  Message sent via/by Savannah
  http://savannah.gnu.org/