/* fft2d.h

Interface to the 'fft2d' DLL.
 
The 'fft2d' DLL exports the functions of the `fft2d library by Takuya Ooura 
<http://www.kurims.kyoto-u.ac.jp/~ooura/fft.html>`_ and implements fast 
3D auto- and cross-correlation functions.

:Author:
  `Christoph Gohlke <http://www.lfd.uci.edu/~gohlke/>`_

:Organization:
  Laboratory for Fluorescence Dynamics, University of California, Irvine

:Version: 2015.07.17

Example
-------
Inplace autocorrelation of a 4x4x4 array::

    #include <stdio.h>
    #include "fft2d.h"
    int main(int argc, char** argv)
    {
        ssize_t strides[] = {4*4*sizeof(float), 4*sizeof(float), sizeof(float)};
        float data[4*4*4];
        for (int i = 0; i < 4*4*4; i++) data[i] = (float)i;
        rfft3d_handle handle = rfft3d_new(4, 4, 4, FFT3D_MODE_TYX | FFT3D_MODE_FCS);
        if (handle == NULL) return -1;
        rfft3d_autocorrelate_ff(handle, data, data, strides);
        rfft3d_del(handle);
        for (int i = 0; i < 4 * 4 * 4;) {
            if (i % 4 == 0) printf("\n");
            printf("%12.9f ", data[i++]);
            if (i % 16 == 0) printf("\n");
        }
        return 0;
    }

License
-------
Copyright (c) 2015, Christoph Gohlke
Copyright (c) 2015, The Regents of the University of California
Produced at the Laboratory for Fluorescence Dynamics.
All rights reserved.

Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:

* Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright
notice, this list of conditions and the following disclaimer in the
documentation and/or other materials provided with the distribution.
* Neither the name of the copyright holders nor the names of any
contributors may be used to endorse or promote products derived
from this software without specific prior written permission.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
POSSIBILITY OF SUCH DAMAGE.
*/

#ifndef FFT2D_H
#define FFT2D_H

#include <stdint.h>

#if defined(_MSC_VER)
#include <BaseTsd.h>
typedef SSIZE_T ssize_t;
#endif

#ifdef FFT2D_USEDLL
# define FFT2D_API extern __declspec(dllimport)
#elif defined FFT2D_EXPORTS
# define FFT2D_API extern __declspec(dllexport)
#else
# define FFT2D_API
#endif

#ifdef USE_FFT2D_WINTHREADS
# ifndef FFT2D_MAX_THREADS
#  define FFT2D_MAX_THREADS 4
# endif
#else
# define FFT2D_MAX_THREADS 1
#endif

#ifdef USE_FFT3D_WINTHREADS
# ifndef FFT3D_MAX_THREADS
#  define FFT3D_MAX_THREADS 4
# endif
#else
# define FFT3D_MAX_THREADS 1
#endif


#define FFT2D_VERSION_STR "2015.07.17"

#define FFT3D_MODE_DEFAULT 0
#define FFT3D_MODE_TYX 1   /* do not center correlation results in time axis */
#define FFT3D_MODE_FCS 2   /* normalize correlation results according to FCS */
#define FFT3D_MODE_CC  4   /* allocate second buffer for cross correlation */

#define FFT2D_ERROR -1
#define FFT2D_VALUE_ERROR -2
#define FFT2D_MEMORY_ERROR -3

/* C++ API */

#ifdef  __cplusplus
    struct rfft3d
    {
    public:
        rfft3d(ssize_t shape0, ssize_t shape1, ssize_t shape2, int mode = 0);
        ~rfft3d();
        void set_mode(int mode) { mode_ = mode; }
        void init_tables() { ip_[0] = 0; }
        template <typename Ti, typename To> void autocorrelate(Ti* data, To* out, ssize_t* strides = NULL);
        template <typename Ti, typename To> void crosscorrelate(Ti* data0, Ti* data1, To* out, ssize_t* strides0 = NULL, ssize_t* strides1 = NULL);
        template <typename Ti, typename To> void fft(Ti* data, To* out, ssize_t* strides = NULL);
        template <typename Ti, typename To> void ifft(Ti* data, To* out, ssize_t* strides = NULL);
    private:
        int mode_;
        ssize_t np_;
        ssize_t n0_, n1_, n2_;  /* shape of input and output arrays */
        double ***a_;  /* input/output data */
        double ***b_;  /* cross correlation input data */
        double *w_;  /* cos/sin table */
        double *t_;  /* work area */
        int *ip_;  /* work area for bit reversal */
        void alloc_();
        void free_();
        template <typename T> double copy_input_(double*** a, T* data, ssize_t* strides);
        template <typename T> void copy_output_(double*** a, T* out, const double scale=1.0, const double offset=0.0);
    };
#else
    struct rfft3d;
#endif

/* C API */

    typedef struct rfft3d* rfft3d_handle;

#ifdef __cplusplus
extern "C" {
#endif

    /* rfft3d class */
    FFT2D_API rfft3d_handle rfft3d_new(ssize_t shape0, ssize_t shape1, ssize_t shape2, int mode);
    FFT2D_API void rfft3d_del(rfft3d_handle);
    FFT2D_API void rfft3d_init_tables(rfft3d_handle);

    FFT2D_API int rfft3d_autocorrelate_dd(rfft3d_handle, double* data, double* out, ssize_t* strides);
    FFT2D_API int rfft3d_autocorrelate_ff(rfft3d_handle, float* data, float* out, ssize_t* strides);
    FFT2D_API int rfft3d_autocorrelate_fd(rfft3d_handle, float* data, double* out, ssize_t* strides);
    FFT2D_API int rfft3d_autocorrelate_hf(rfft3d_handle, int16_t* data, float* out, ssize_t* strides);
    FFT2D_API int rfft3d_autocorrelate_hd(rfft3d_handle, int16_t* data, double* out, ssize_t* strides);

    FFT2D_API int rfft3d_crosscorrelate_dd(rfft3d_handle, double* data0, double* data1, double* out, ssize_t* strides0, ssize_t* strides1);
    FFT2D_API int rfft3d_crosscorrelate_ff(rfft3d_handle, float* data0, float* data1, float* out, ssize_t* strides0, ssize_t* strides1);
    FFT2D_API int rfft3d_crosscorrelate_fd(rfft3d_handle, float* data0, float* data1, double* out, ssize_t* strides0, ssize_t* strides1);
    FFT2D_API int rfft3d_crosscorrelate_hf(rfft3d_handle, int16_t* data0, int16_t* data1, float* out, ssize_t* strides0, ssize_t* strides1);
    FFT2D_API int rfft3d_crosscorrelate_hd(rfft3d_handle, int16_t* data0, int16_t* data1, double* out, ssize_t* strides0, ssize_t* strides1);

    /* FFT2D library functions */
    FFT2D_API int *alloc_1d_int(int n1);
    FFT2D_API void free_1d_int(int *i);
    FFT2D_API double *alloc_1d_double(int n1);
    FFT2D_API void free_1d_double(double *d);
    FFT2D_API int **alloc_2d_int(int n1, int n2);
    FFT2D_API void free_2d_int(int **ii);
    FFT2D_API double **alloc_2d_double(int n1, int n2);
    FFT2D_API void free_2d_double(double **dd);
    FFT2D_API int ***alloc_3d_int(int n1, int n2, int n3);
    FFT2D_API void free_3d_int(int ***iii);
    FFT2D_API double ***alloc_3d_double(int n1, int n2, int n3);
    FFT2D_API void free_3d_double(double ***ddd);

    FFT2D_API void cdft(int, int, double *, int *, double *);
    FFT2D_API void rdft(int, int, double *, int *, double *);
    FFT2D_API void ddct(int, int, double *, int *, double *);
    FFT2D_API void ddst(int, int, double *, int *, double *);
    FFT2D_API void dfct(int, double *, double *, int *, double *);
    FFT2D_API void dfst(int, double *, double *, int *, double *);

    FFT2D_API void cdft2d(int, int, int, double **, double *, int *, double *);
    FFT2D_API void rdft2d(int, int, int, double **, double *, int *, double *);
    FFT2D_API void rdft2dsort(int, int, int, double **);
    FFT2D_API void ddct2d(int, int, int, double **, double *, int *, double *);
    FFT2D_API void ddst2d(int, int, int, double **, double *, int *, double *);

    FFT2D_API void cdft3d(int, int, int, int, double ***, double *, int *, double *);
    FFT2D_API void rdft3d(int, int, int, int, double ***, double *, int *, double *);
    FFT2D_API void rdft3dsort(int, int, int, int, double ***);
    FFT2D_API void ddct3d(int, int, int, int, double ***, double *, int *, double *);
    FFT2D_API void ddst3d(int, int, int, int, double ***, double *, int *, double *);

    FFT2D_API void ddct8x8s(int isgn, double **a);
    FFT2D_API void ddct16x16s(int isgn, double **a);

#ifdef __cplusplus
}
#endif

#endif /* FFT2D_H */