#include <Python.h>
#include <numpy/arrayobject.h>

#define NBDIM 3

/*  wrapped refine function */
static PyObject* refine(PyObject* self, PyObject* args)
{
	PyArrayObject* in_array;
	int rx, ry, rz;
	PyObject* out_array;
	NpyIter* in_iter;
	NpyIter* out_iter;
	NpyIter_IterNextFunc* in_iternext;
	NpyIter_IterNextFunc* out_iternext;

	/*  parse single numpy array argument */
	if (!PyArg_ParseTuple(args, "O!iii", &PyArray_Type, &in_array, &rz, &ry, &rx))
		return NULL;
	if (PyArray_NDIM(in_array) != NBDIM)
		return NULL;

	PyArray_Descr out_descr;
	out_descr=*PyArray_DESCR(in_array);



	/*  construct the output array, like the input array */
	out_array = PyArray_NewFromDescr(&PyArray_Type, out_descr, NBDIM, dims, NULL, data, NPY_ARRAY_OWNDATA, NULL);
	if (out_array == NULL)
		return NULL;

	/*  create the iterators */
	in_iter = NpyIter_New(in_array, NPY_ITER_READONLY, NPY_KEEPORDER, NPY_NO_CASTING, NULL);
	if (in_iter == NULL)
		goto fail;

	out_iter = NpyIter_New((PyArrayObject*)out_array, NPY_ITER_READWRITE, NPY_KEEPORDER, NPY_NO_CASTING, NULL);
	if (out_iter == NULL) {
		NpyIter_Deallocate(in_iter);
		goto fail;
	}

	in_iternext = NpyIter_GetIterNext(in_iter, NULL);
	out_iternext = NpyIter_GetIterNext(out_iter, NULL);
	if (in_iternext == NULL || out_iternext == NULL) {
		NpyIter_Deallocate(in_iter);
		NpyIter_Deallocate(out_iter);
		goto fail;
	}



    double** in_dataptr = (double**) NpyIter_GetDataPtrArray(in_iter);
    double** out_dataptr = (double**) NpyIter_GetDataPtrArray(out_iter);

    /*  iterate over the arrays */
    do {
        **out_dataptr = cos(**in_dataptr);
    } while(in_iternext(in_iter) && out_iternext(out_iter));

    /*  clean up and return the result */
    NpyIter_Deallocate(in_iter);
    NpyIter_Deallocate(out_iter);
    Py_INCREF(out_array);
    return out_array;

    /*  in case bad things happen */
    fail:
        Py_XDECREF(out_array);
        return NULL;
}

/*  define functions in module */
static PyMethodDef RefineMethods[] =
{
     {"refine", refine, METH_VARARGS, "Refines a numpy array"},
     {NULL, NULL, 0, NULL}
};

static struct PyModuleDef crefineDef =
{
    PyModuleDef_HEAD_INIT,
    "refine", "",
    -1,
    RefineMethods
};


import_array();
PyMODINIT_FUNC
PyInit_refine(void)
{
	return PyModule_Create(&crefineDef);
}