Merge pull request #178 from sdss/fixsetup

Proper installation of fast median extension

README.md

@@ -63,18 +63,7 @@ To install the DRP along with its dependencies, you need to run the following st
     cd lvmdrp
     ```
 
-3. Optional, but recommended: compile c-code for better performance
-
-    ```bash
-    pushd lvmdrp/python/cextern/fast_median/src
-    make
-    popd
-    ```
-    This should have created a fast_median.so (Linux) or fast_median.dylib (MacOS) in the src directory.
-    Set environment variable `LVMDRP_LIB_DIR` if you move the shared library. If you leave it in src/, no need to set this.
-    We will automate this step eventually ;-)
-
-4. Install the DRP package in the current python environment (see [contributing](#contributing-to-lvm-drp-development) section below for a replacement of this step). Make sure you're back in the lvmdrp directory.
+3. Install the DRP package in the current python environment (see [contributing](#contributing-to-lvm-drp-development) section below for a replacement of this step). Make sure you're back in the lvmdrp directory.
 
     ```bash
     cd lvmdrp

cextern/README.txt

@@ -0,0 +1,7 @@
+C++ code to generate fast_median extension goes here.
+
+fast_median.cpp, fast_median.hpp: a very fast median filter for 1d and 2d ndarrays
+with float or double data.
+python/lvmdrp/external/fast_median.py: python interface for fast_median extension
+
+

python/cextern/fast_median/src/fast_median.cpp → cextern/fast_median/src/fast_median.cpp

@@ -5,12 +5,14 @@
 #include <stdexcept>
 #include <stdint.h>
 
+
 #ifndef __arm64__
 #include <x86intrin.h>
 #endif
 
 // see https://github.com/suomela/mf2d
 
+#include <Python.h>
 #include "fast_median.hpp"
 
 const uint64_t ONE64 = 1;
@@ -667,3 +669,115 @@ void median_filter_1d_double(int nx, int box_x, int blockhint, const double *in,
     return median_filter_1d<double>(nx, hx, blockhint, in, out);
 }
 
+
+// Wrapper for median_filter_1d_float
+static PyObject* py_median_filter_1d_float(PyObject* self, PyObject* args) {
+    Py_buffer input_buf, output_buf;
+    int nx, box_x, blockhint;
+
+    if (!PyArg_ParseTuple(args, "iiiy*y*", &nx, &box_x, &blockhint, &input_buf, &output_buf)) {
+        PyErr_SetString(PyExc_ValueError, "Failed to parse arguments. Expected: int, int, int, buffer, buffer.");
+        return NULL;
+    }
+
+    if (input_buf.len != output_buf.len || nx * sizeof(float) != input_buf.len) {
+        PyErr_SetString(PyExc_ValueError, "Buffer sizes do not match input dimensions.");
+        return NULL;
+    }
+
+    median_filter_1d_float(nx, box_x, blockhint, (const float*)input_buf.buf, (float*)output_buf.buf);
+
+    PyBuffer_Release(&input_buf);
+    PyBuffer_Release(&output_buf);
+
+    Py_RETURN_NONE;
+}
+
+// Wrapper for median_filter_1d_double
+static PyObject* py_median_filter_1d_double(PyObject* self, PyObject* args) {
+    Py_buffer input_buf, output_buf;
+    int nx, box_x, blockhint;
+
+    if (!PyArg_ParseTuple(args, "iiiy*y*", &nx, &box_x, &blockhint, &input_buf, &output_buf)) {
+        PyErr_SetString(PyExc_ValueError, "Failed to parse arguments. Expected: int, int, int, buffer, buffer.");
+        return NULL;
+    }
+
+    if (input_buf.len != output_buf.len || nx * sizeof(double) != input_buf.len) {
+        PyErr_SetString(PyExc_ValueError, "Buffer sizes do not match input dimensions.");
+        return NULL;
+    }
+
+    median_filter_1d_double(nx, box_x, blockhint, (const double*)input_buf.buf, (double*)output_buf.buf);
+
+    PyBuffer_Release(&input_buf);
+    PyBuffer_Release(&output_buf);
+
+    Py_RETURN_NONE;
+}
+
+// Wrapper for median_filter_2d_float
+static PyObject* py_median_filter_2d_float(PyObject* self, PyObject* args) {
+    Py_buffer input_buf, output_buf;
+    int nx, ny, box_x, box_y, blockhint;
+
+    if (!PyArg_ParseTuple(args, "iiiiiy*y*", &nx, &ny, &box_x, &box_y, &blockhint, &input_buf, &output_buf)) {
+        PyErr_SetString(PyExc_ValueError, "Failed to parse arguments. Expected: int, int, int, int, int, buffer, buffer.");
+        return NULL;
+    }
+
+    if (input_buf.len != output_buf.len || nx * ny * sizeof(float) != input_buf.len) {
+        PyErr_SetString(PyExc_ValueError, "Buffer sizes do not match input dimensions.");
+        return NULL;
+    }
+
+    median_filter_2d_float(nx, ny, box_x, box_y, blockhint, (const float*)input_buf.buf, (float*)output_buf.buf);
+
+    PyBuffer_Release(&input_buf);
+    PyBuffer_Release(&output_buf);
+
+    Py_RETURN_NONE;
+}
+
+// Wrapper for median_filter_2d_double
+static PyObject* py_median_filter_2d_double(PyObject* self, PyObject* args) {
+    Py_buffer input_buf, output_buf;
+    int nx, ny, box_x, box_y, blockhint;
+
+    if (!PyArg_ParseTuple(args, "iiiiiy*y*", &nx, &ny, &box_x, &box_y, &blockhint, &input_buf, &output_buf)) {
+        PyErr_SetString(PyExc_ValueError, "Failed to parse arguments. Expected: int, int, int, int, int, buffer, buffer.");
+        return NULL;
+    }
+
+    if (input_buf.len != output_buf.len || nx * ny * sizeof(double) != input_buf.len) {
+        PyErr_SetString(PyExc_ValueError, "Buffer sizes do not match input dimensions.");
+        return NULL;
+    }
+
+    median_filter_2d_double(nx, ny, box_x, box_y, blockhint, (const double*)input_buf.buf, (double*)output_buf.buf);
+
+    PyBuffer_Release(&input_buf);
+    PyBuffer_Release(&output_buf);
+
+    Py_RETURN_NONE;
+}
+
+static PyMethodDef FastMedianMethods[] = {
+    {"median_filter_1d_float", py_median_filter_1d_float, METH_VARARGS, "1D median filter for float arrays"},
+    {"median_filter_1d_double", py_median_filter_1d_double, METH_VARARGS, "1D median filter for double arrays"},
+    {"median_filter_2d_float", py_median_filter_2d_float, METH_VARARGS, "2D median filter for float arrays"},
+    {"median_filter_2d_double", py_median_filter_2d_double, METH_VARARGS, "2D median filter for double arrays"},
+    {NULL, NULL, 0, NULL}  // Sentinel
+};
+
+static struct PyModuleDef fast_median_module = {
+    PyModuleDef_HEAD_INIT,
+    "fast_median",                  // Module name
+    "Fast median filters",          // Module docstring
+    -1,                             // Size of per-interpreter state or -1 for global
+    FastMedianMethods               // Method table
+};
+
+extern "C" PyMODINIT_FUNC PyInit_fast_median(void) {
+    return PyModule_Create(&fast_median_module);
+}

python/lvmdrp/core/image.py

@@ -25,7 +25,7 @@
 from lvmdrp.core.tracemask import TraceMask
 from lvmdrp.core.spectrum1d import Spectrum1D, _normalize_peaks, _fiber_cc_match, _cross_match, _spec_from_lines, align_blocks
 
-from cextern.fast_median.fast_median import fast_median_filter_2d
+from lvmdrp.external.fast_median import fast_median_filter_2d
 
 from lvmdrp import __version__ as drpver