ENH: enable complex dtype support unconditionally (#790)

* WIP: ENH: enable complex dtype support unconditionally

This also gets rid of usage of the deprecated IF macro in Cython code,
which resolves a bunch of build warnings.

* BLD/ENH: support complex dtypes on MSVC, mirroring NumPy's approach

MSVC has no C99 `_Complex`; it provides the struct-based
`_Fcomplex`/`_Dcomplex` types, which have no arithmetic operators. Handle
this the way NumPy does (see npy_common.h/npy_math.h and the complex
check in numpy/_core/meson.build):

- common.h: typedef switch (`_Fcomplex`/`_Dcomplex` on MSVC, C99 complex
  elsewhere) plus inline accessors that read via `creal` and write via
  pointer punning, valid for every representation. A
  PYWT_TEST_STRUCT_COMPLEX define forces NumPy's struct branch so the
  MSVC code path can be compiled and tested with any compiler.
- convolution.template.c (the only complex arithmetic in the C core):
  replace operators with type-generic helpers - ZERO/ADD/SUB/MUL and
  ADD_PRODUCT(s, f, x) for `s += f * x`. The multiplier is always real,
  so real and imaginary parts stay independent and codegen on
  C99-complex compilers is unchanged.
- Cython boundary: declare float_complex/double_complex as opaque extern
  typedefs in c_wt.pxd and cast explicitly at call sites, instead of
  forcing -DCYTHON_CCOMPLEX=1; Cython >=3.2.5 picks the right complex
  representation per compiler.
- meson: replace the have_c99_complex probe with NumPy's complex.h/type
  checks, and drop the _c99_config machinery plus the remaining
  real/imag-splitting fallbacks in _swt.py/_multidim.py.

The full test suite passes both normally and with the struct
representation forced (emulating MSVC semantics on Linux).

* CI: test the MSVC-style struct complex path on linux-aarch64

Build with -DPYWT_TEST_STRUCT_COMPLEX so the operator-free complex code
path that MSVC relies on is compiled and tested on every PR, not only in
the Windows wheel builds.

* DOC: clarify the PYWT_TEST_STRUCT_COMPLEX comment in common.h

* BLD: document why the Cython copy-tree must stay in meson.build

The generated config.pxi is gone, but the copies still provide the
rebuild dependency on .pxd files: Cython's -M depfiles (which Meson's
native Cython support relies on) do not record module-namespace cimports
such as 'from . cimport c_wt', so dropping the copy-tree would cause
stale builds after .pxd edits.

Author: rgommers

.github/workflows/tests.yml

@@ -41,6 +41,12 @@ jobs:
           # Linux arm64
           - runs-on: ubuntu-22.04-arm
             python-version: "3.14"
+          # Linux arm64, MSVC-style struct complex types forced, to exercise
+          # the operator-free complex code path on a C99-complex compiler
+          # (see PYWT_TEST_STRUCT_COMPLEX in pywt/_extensions/c/common.h)
+          - runs-on: ubuntu-22.04-arm
+            python-version: "3.14"
+            OPTIONS_NAME: "struct-complex"
           # Linux amd64
           - runs-on: ubuntu-latest
             python-version: "3.12"
@@ -100,6 +106,8 @@ jobs:
             elif [ "${REFGUIDE_CHECK}" == "1" ]; then
                 pip install sphinx numpydoc scipy-doctest
                 pip install . -v
+            elif [ "${OPTIONS_NAME}" == "struct-complex" ]; then
+                pip install . -v -Csetup-args=-Dc_args=-DPYWT_TEST_STRUCT_COMPLEX
             else
                 pip install . -v
             fi

.gitignore

@@ -36,8 +36,6 @@ doc/source/regression/*.ipynb
 # Project working files
 # Expanded Cython
 pywt/_extensions/*.[ch]
-pywt/_extensions/_c99_config.py
-pywt/_extensions/config.pxi
 cythonize.dat
 pywt/version.py
 build.log

meson.build

@@ -11,7 +11,7 @@ project(
   ],
 )
 
-cython = find_program('cython')
+cython = find_program('cython', version: '>=3.2.5')
 
 py = import('python').find_installation(pure: false)
 py_dep = py.dependency()

pyproject.toml

@@ -8,7 +8,7 @@
 build-backend = "mesonpy"
 requires = [
     "meson-python>=0.18.0",
-    "Cython>=3.1.3",
+    "Cython>=3.2.5",
 
     # numpy requirement for wheel builds for distribution on PyPI - building
     # against 2.x yields wheels that are also compatible with numpy 1.x at

pywt/_c99_config.py.in

@@ -2,7 +2,6 @@
 
 import numpy as np
 
-from ._c99_config import _have_c99_complex
 from ._extensions._dwt import downcoef as _downcoef
 from ._extensions._dwt import dwt_axis, dwt_single, idwt_axis, idwt_single
 from ._extensions._dwt import dwt_coeff_len as _dwt_coeff_len
@@ -161,12 +160,6 @@ def dwt(data, wavelet, mode='symmetric', axis=-1):
     array([-0.70710678, -0.70710678, -0.70710678])
 
     """
-    if not _have_c99_complex and np.iscomplexobj(data):
-        data = np.asarray(data)
-        cA_r, cD_r = dwt(data.real, wavelet, mode, axis)
-        cA_i, cD_i = dwt(data.imag, wavelet, mode, axis)
-        return (cA_r + 1j*cA_i, cD_r + 1j*cD_i)
-
     # accept array_like input; make a copy to ensure a contiguous array
     dt = _check_dtype(data)
     data = np.asarray(data, dtype=dt, order='C')
@@ -241,17 +234,6 @@ def idwt(cA, cD, wavelet, mode='symmetric', axis=-1):
         raise ValueError("At least one coefficient parameter must be "
                          "specified.")
 
-    # for complex inputs: compute real and imaginary separately then combine
-    if not _have_c99_complex and (np.iscomplexobj(cA) or np.iscomplexobj(cD)):
-        if cA is None:
-            cD = np.asarray(cD)
-            cA = np.zeros_like(cD)
-        elif cD is None:
-            cA = np.asarray(cA)
-            cD = np.zeros_like(cA)
-        return (idwt(cA.real, cD.real, wavelet, mode, axis) +
-                1j*idwt(cA.imag, cD.imag, wavelet, mode, axis))
-
     if cA is not None:
         dt = _check_dtype(cA)
         cA = np.asarray(cA, dtype=dt, order='C')
@@ -328,9 +310,6 @@ def downcoef(part, data, wavelet, mode='symmetric', level=1):
     upcoef
 
     """
-    if not _have_c99_complex and np.iscomplexobj(data):
-        return (downcoef(part, data.real, wavelet, mode, level) +
-                1j*downcoef(part, data.imag, wavelet, mode, level))
     # accept array_like input; make a copy to ensure a contiguous array
     dt = _check_dtype(data)
     data = np.asarray(data, dtype=dt, order='C')
@@ -387,9 +366,6 @@ def upcoef(part, coeffs, wavelet, level=1, take=0):
     array([ 1.,  2.,  3.,  4.,  5.,  6.])
 
     """
-    if not _have_c99_complex and np.iscomplexobj(coeffs):
-        return (upcoef(part, coeffs.real, wavelet, level, take) +
-                1j*upcoef(part, coeffs.imag, wavelet, level, take))
     # accept array_like input; make a copy to ensure a contiguous array
     dt = _check_dtype(coeffs)
     coeffs = np.asarray(coeffs, dtype=dt, order='C')