# HG changeset patch
# User Markus Meisinger <chloros2@gmx.de>
# Date 1640431561 -3600
#      Sat Dec 25 12:26:01 2021 +0100
# Node ID 5ae7b66c4959a74bf976b5e3c581ad72c46e40bb
# Parent  b65d423b12a6a391867663ec6b51c27c8a6f99bb
ifft, ifft2, ifftn: Return the real part if the symflag flag is specified (bug #43742).

* fft.cc, fft2.cc, fftn.cc: Peel "symmetric" keyword and return real part of ifft().

diff -r b65d423b12a6 -r 5ae7b66c4959 libinterp/corefcn/fft.cc
--- a/libinterp/corefcn/fft.cc	Sun Dec 19 20:42:45 2021 +0100
+++ b/libinterp/corefcn/fft.cc	Sat Dec 25 12:26:01 2021 +0100
@@ -41,8 +41,9 @@
 do_fft (const octave_value_list& args, const char *fcn, int type)
 {
   int nargin = args.length ();
+  int numeric_argin = nargin;
 
-  if (nargin < 1 || nargin > 3)
+  if (nargin < 1 || nargin > 4)
     print_usage ();
 
   octave_value retval;
@@ -51,8 +52,26 @@
   dim_vector dims = arg.dims ();
   int ndims = dims.ndims ();
   int dim = -1;
+  bool have_symmetric = false;
 
-  if (nargin > 1)
+  // Peel "symflag"
+  if ((nargin > 1) && (! args(nargin - 1).isempty ()) &&
+      (! args(nargin - 1).isnumeric ()))
+    {
+      std::string strval = args(nargin - 1).xstring_value ("must be a string");
+      if (strval == "symmetric")
+        {
+          have_symmetric = true;
+        }
+      else if (strval != "nonsymmetric")
+        {
+          error ("%s: unknown symflag", fcn);
+        }
+      numeric_argin = nargin - 1;
+    }
+
+  // Get numeric arguments
+  if (numeric_argin > 1)
     {
       if (! args(1).isempty ())
         {
@@ -66,7 +85,7 @@
         }
     }
 
-  if (nargin > 2)
+  if (numeric_argin > 2)
     {
       double dval = args(2).double_value ();
       if (math::isnan (dval))
@@ -122,34 +141,54 @@
     {
       if (arg.isreal ())
         {
-          FloatNDArray nda = arg.float_array_value ();
+          FloatNDArray in = arg.float_array_value ();
+
+          in.resize (dims, 0.0);
+          FloatComplexNDArray out = (type != 0 ? in.ifourier (dim) : in.fourier (dim));
 
-          nda.resize (dims, 0.0);
-          retval = (type != 0 ? nda.ifourier (dim) : nda.fourier (dim));
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
       else
         {
-          FloatComplexNDArray cnda = arg.float_complex_array_value ();
+          FloatComplexNDArray in = arg.float_complex_array_value ();
+
+          in.resize (dims, 0.0);
+          FloatComplexNDArray out = (type != 0 ? in.ifourier (dim) : in.fourier (dim));
 
-          cnda.resize (dims, 0.0);
-          retval = (type != 0 ? cnda.ifourier (dim) : cnda.fourier (dim));
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
     }
   else
     {
       if (arg.isreal ())
         {
-          NDArray nda = arg.array_value ();
+          NDArray in = arg.array_value ();
+
+          in.resize (dims, 0.0);
+          ComplexNDArray out = (type != 0 ? in.ifourier (dim) : in.fourier (dim));
 
-          nda.resize (dims, 0.0);
-          retval = (type != 0 ? nda.ifourier (dim) : nda.fourier (dim));
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
       else if (arg.iscomplex ())
         {
-          ComplexNDArray cnda = arg.complex_array_value ();
+          ComplexNDArray in = arg.complex_array_value ();
+
+          in.resize (dims, 0.0);
+          ComplexNDArray out = (type != 0 ? in.ifourier (dim) : in.fourier (dim));
 
-          cnda.resize (dims, 0.0);
-          retval = (type != 0 ? cnda.ifourier (dim) : cnda.fourier (dim));
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
       else
         err_wrong_type_arg (fcn, arg);
@@ -227,6 +266,7 @@
 @deftypefn  {} {} ifft (@var{x})
 @deftypefnx {} {} ifft (@var{x}, @var{n})
 @deftypefnx {} {} ifft (@var{x}, @var{n}, @var{dim})
+@deftypefnx {} {} ifft (@dots{}, @var{symflag})
 Compute the inverse discrete Fourier transform of @var{x}
 using a Fast Fourier Transform (FFT) algorithm.
 
@@ -244,6 +284,10 @@
 
 If called with three arguments, @var{dim} is an integer specifying the
 dimension of the matrix along which the inverse FFT is performed.
+
+If the FT satisfies Hermitian redundancy, it's inverse is real-valued. In
+order to force a real-valued output the keyword @code{@var{symmetric}} can
+be passed as an argument.
 @seealso{fft, ifft2, ifftn, fftw}
 @end deftypefn */)
 {
diff -r b65d423b12a6 -r 5ae7b66c4959 libinterp/corefcn/fft2.cc
--- a/libinterp/corefcn/fft2.cc	Sun Dec 19 20:42:45 2021 +0100
+++ b/libinterp/corefcn/fft2.cc	Sat Dec 25 12:26:01 2021 +0100
@@ -43,16 +43,35 @@
 do_fft2 (const octave_value_list& args, const char *fcn, int type)
 {
   int nargin = args.length ();
+  int numeric_argin = nargin;
 
-  if (nargin < 1 || nargin > 3)
+  if (nargin < 1 || nargin > 4)
     print_usage ();
 
   octave_value retval;
   octave_value arg = args(0);
   dim_vector dims = arg.dims ();
   octave_idx_type n_rows = -1;
+  bool have_symmetric = false;
 
-  if (nargin > 1)
+  // Peel "symflag"
+  if ((nargin > 1) && (! args(nargin - 1).isempty ()) &&
+      (! args(nargin - 1).isnumeric ()))
+    {
+      std::string strval = args(nargin - 1).xstring_value ("must be a string");
+      if (strval == "symmetric")
+        {
+          have_symmetric = true;
+        }
+      else if (strval != "nonsymmetric")
+        {
+          error ("%s: unknown symflag", fcn);
+        }
+      numeric_argin = nargin - 1;
+    }
+
+  // Get numeric arguments
+  if (numeric_argin > 1)
     {
       double dval = args(1).double_value ();
       if (math::isnan (dval))
@@ -64,7 +83,7 @@
     }
 
   octave_idx_type n_cols = -1;
-  if (nargin > 2)
+  if (numeric_argin > 2)
     {
       double dval = args(2).double_value ();
       if (math::isnan (dval))
@@ -101,34 +120,54 @@
     {
       if (arg.isreal ())
         {
-          FloatNDArray nda = arg.float_array_value ();
+          FloatNDArray in = arg.float_array_value ();
+
+          in.resize (dims, 0.0);
+          FloatComplexNDArray out = (type != 0 ? in.ifourier2d () : in.fourier2d ());
 
-          nda.resize (dims, 0.0);
-          retval = (type != 0 ? nda.ifourier2d () : nda.fourier2d ());
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
       else
         {
-          FloatComplexNDArray cnda = arg.float_complex_array_value ();
+          FloatComplexNDArray in = arg.float_complex_array_value ();
+
+          in.resize (dims, 0.0);
+          FloatComplexNDArray out = (type != 0 ? in.ifourier2d () : in.fourier2d ());
 
-          cnda.resize (dims, 0.0);
-          retval = (type != 0 ? cnda.ifourier2d () : cnda.fourier2d ());
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
     }
   else
     {
       if (arg.isreal ())
         {
-          NDArray nda = arg.array_value ();
+          NDArray in = arg.array_value ();
+
+          in.resize (dims, 0.0);
+          ComplexNDArray out = (type != 0 ? in.ifourier2d () : in.fourier2d ());
 
-          nda.resize (dims, 0.0);
-          retval = (type != 0 ? nda.ifourier2d () : nda.fourier2d ());
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
       else if (arg.iscomplex ())
         {
-          ComplexNDArray cnda = arg.complex_array_value ();
+          ComplexNDArray in = arg.complex_array_value ();
+
+          in.resize (dims, 0.0);
+          ComplexNDArray out = (type != 0 ? in.ifourier2d () : in.fourier2d ());
 
-          cnda.resize (dims, 0.0);
-          retval = (type != 0 ? cnda.ifourier2d () : cnda.fourier2d ());
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
       else
         err_wrong_type_arg (fcn, arg);
@@ -161,6 +200,7 @@
        doc: /* -*- texinfo -*-
 @deftypefn  {} {} ifft2 (@var{A})
 @deftypefnx {} {} ifft2 (@var{A}, @var{m}, @var{n})
+@deftypefnx {} {} ifft2 (@dots{}, @var{symflag})
 Compute the inverse two-dimensional discrete Fourier transform of @var{A}
 using a Fast Fourier Transform (FFT) algorithm.
 
@@ -170,6 +210,10 @@
 
 If @var{A} is a multi-dimensional matrix, each two-dimensional sub-matrix
 of @var{A} is treated separately.
+
+If the FT satisfies Hermitian redundancy, it's inverse is real-valued. In
+order to force a real-valued output the keyword @code{@var{symmetric}} can
+be passed as an argument.
 @seealso{fft2, ifft, ifftn, fftw}
 @end deftypefn */)
 {
diff -r b65d423b12a6 -r 5ae7b66c4959 libinterp/corefcn/fftn.cc
--- a/libinterp/corefcn/fftn.cc	Sun Dec 19 20:42:45 2021 +0100
+++ b/libinterp/corefcn/fftn.cc	Sat Dec 25 12:26:01 2021 +0100
@@ -43,19 +43,38 @@
 do_fftn (const octave_value_list& args, const char *fcn, int type)
 {
   int nargin = args.length ();
+  int numeric_argin = nargin;
 
-  if (nargin < 1 || nargin > 2)
+  if (nargin < 1 || nargin > 3)
     print_usage ();
 
   octave_value retval;
   octave_value arg = args(0);
   dim_vector dims = arg.dims ();
+  bool have_symmetric = false;
 
   for (int i = 0; i < dims.ndims (); i++)
     if (dims(i) < 0)
       return retval;
 
-  if (nargin > 1)
+  // Peel "symflag"
+  if ((nargin > 1) && (! args(nargin - 1).isempty ()) &&
+      (! args(nargin - 1).isnumeric ()))
+    {
+      std::string strval = args(nargin - 1).xstring_value ("must be a string");
+      if (strval == "symmetric")
+        {
+          have_symmetric = true;
+        }
+      else if (strval != "nonsymmetric")
+        {
+          error ("%s: unknown symflag", fcn);
+        }
+      numeric_argin = nargin - 1;
+    }
+
+  // Get numeric arguments
+  if (numeric_argin > 1)
     {
       Matrix val = args(1).xmatrix_value ("%s: SIZE must be a vector of length dim", fcn);
 
@@ -88,34 +107,54 @@
     {
       if (arg.isreal ())
         {
-          FloatNDArray nda = arg.float_array_value ();
+          FloatNDArray in = arg.float_array_value ();
+
+          in.resize (dims, 0.0);
+          FloatComplexNDArray out = (type != 0 ? in.ifourierNd () : in.fourierNd ());
 
-          nda.resize (dims, 0.0);
-          retval = (type != 0 ? nda.ifourierNd () : nda.fourierNd ());
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
       else
         {
-          FloatComplexNDArray cnda = arg.float_complex_array_value ();
+          FloatComplexNDArray in = arg.float_complex_array_value ();
+
+          in.resize (dims, 0.0);
+          FloatComplexNDArray out = (type != 0 ? in.ifourierNd () : in.fourierNd ());
 
-          cnda.resize (dims, 0.0);
-          retval = (type != 0 ? cnda.ifourierNd () : cnda.fourierNd ());
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
     }
   else
     {
       if (arg.isreal ())
         {
-          NDArray nda = arg.array_value ();
+          NDArray in = arg.array_value ();
+
+          in.resize (dims, 0.0);
+          ComplexNDArray out = (type != 0 ? in.ifourierNd () : in.fourierNd ());
 
-          nda.resize (dims, 0.0);
-          retval = (type != 0 ? nda.ifourierNd () : nda.fourierNd ());
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
       else if (arg.iscomplex ())
         {
-          ComplexNDArray cnda = arg.complex_array_value ();
+          ComplexNDArray in = arg.complex_array_value ();
+
+          in.resize (dims, 0.0);
+          ComplexNDArray out = (type != 0 ? in.ifourierNd () : in.fourierNd ());
 
-          cnda.resize (dims, 0.0);
-          retval = (type != 0 ? cnda.ifourierNd () : cnda.fourierNd ());
+          if (type && have_symmetric)
+            retval = real (out);
+          else
+            retval = out;
         }
       else
         err_wrong_type_arg (fcn, arg);
@@ -147,6 +186,7 @@
        doc: /* -*- texinfo -*-
 @deftypefn  {} {} ifftn (@var{A})
 @deftypefnx {} {} ifftn (@var{A}, @var{size})
+@deftypefnx {} {} ifftn (@dots{}, @var{symflag})
 Compute the inverse N-dimensional discrete Fourier transform of @var{A}
 using a Fast Fourier Transform (FFT) algorithm.
 
@@ -156,6 +196,10 @@
 truncated prior to performing the inverse FFT@.  Otherwise, if an element of
 @var{size} is larger than the corresponding dimension then @var{A} is
 resized and padded with zeros.
+
+If the FT satisfies Hermitian redundancy, it's inverse is real-valued. In
+order to force a real-valued output the keyword @code{@var{symmetric}} can
+be passed as an argument.
 @seealso{fftn, ifft, ifft2, fftw}
 @end deftypefn */)
 {
