JP-3698: Fixed Memory Leak in Ramp Fitting C-Extension (#281)

CHANGES.rst

@@ -16,6 +16,11 @@ Bug Fixes
 
 - 
 
+ramp_fitting
+~~~~~~~~~~~~
+
+- Fixed memory leak in C-extension.[#281]
+
 1.8.0 (2024-08-14)
 ==================
 

src/stcal/ramp_fitting/src/slope_fitter.c

@@ -3,9 +3,13 @@
 
 #include <stdlib.h>
 #include <math.h>
-#include <stdio.h>
 #include <stdbool.h>
+#include <stdio.h>
+#include <stdlib.h>
 #include <string.h>
+#include <sys/resource.h>
+#include <time.h>
+#include <unistd.h>
 
 #include <numpy/arrayobject.h>
 #include <numpy/npy_math.h>
@@ -705,6 +709,30 @@ is_pix_in_list(struct pixel_ramp * pr)
  return 0;
 }
 
+static inline long long
+print_pid_info(long long prev, int line, char * label) {
+ struct rusage res_usage;
+ long long now_time = (long long)time(NULL);
+ long long mem_usage = -1;
+ long long diff = 0;
+ pid_t pid = getpid();
+ // dbg_ols_print("PID: %d\n", pid);
+
+ getrusage(RUSAGE_SELF, &res_usage);
+ mem_usage = res_usage.ru_maxrss;
+ if (prev > 0) {
+ diff = mem_usage - prev;
+ dbg_ols_print("[%d] time: %lld, Mem: %lld, diff: %lld, prev: %lld, pid: %d '%s'\n",
+ line, now_time, mem_usage, diff, prev, pid, label);
+ } else {
+ dbg_ols_print("[%d] time: %lld, Mem: %lld, diff: %lld, pid: %d '%s'\n",
+ line, now_time, mem_usage, diff, pid, label);
+ }
+
+ return mem_usage;
+}
+
+
 /* ------------------------------------------------------------------------- */
 /* Module Functions */
 /* ------------------------------------------------------------------------- */
@@ -899,23 +927,29 @@ clean_ramp_data(
  struct simple_ll_node * current;
  struct simple_ll_node * next;
 
- if (NULL == rd->segs) {
- return; /* Nothing to do. */
- }
-
- /*
- * For each pixel, check to see if there is any allocated
- * memory for the linked list of ramp segments and free them.
- */
- for (idx=0; idx < rd->cube_sz; ++idx) {
- current = rd->segs[idx];
- if (current) {
- next = current->flink;
- SET_FREE(current);
- current = next;
+ Py_XDECREF(rd->data);
+ Py_XDECREF(rd->groupdq);
+ Py_XDECREF(rd->err);
+ Py_XDECREF(rd->pixeldq);
+ Py_XDECREF(rd->zframe);
+ Py_XDECREF(rd->dcurrent);
+
+ if (rd->segs) {
+ /*
+ * For each pixel, check to see if there is any allocated
+ * memory for the linked list of ramp segments and free them.
+ */
+ for (idx=0; idx < rd->cube_sz; ++idx) {
+ current = rd->segs[idx];
+ if (current) {
+ next = current->flink;
+ SET_FREE(current);
+ current = next;
+ }
  }
  }
  SET_FREE(rd->segs);
+ SET_FREE(rd->pedestal);
 }
 
 /*
@@ -1134,6 +1168,9 @@ create_opt_res(
  return 0;
 
 FAILED_ALLOC:
+ PyErr_SetString(PyExc_MemoryError, msg);
+ err_ols_print("%s\n", msg);
+
  Py_XDECREF(opt_res->slope);
  Py_XDECREF(opt_res->sigslope);
  Py_XDECREF(opt_res->var_p);
@@ -1240,6 +1277,9 @@ create_rate_product(
  return 0;
 
 FAILED_ALLOC:
+ PyErr_SetString(PyExc_MemoryError, msg);
+ err_ols_print("%s\n", msg);
+
  Py_XDECREF(rate->slope);
  Py_XDECREF(rate->dq);
  Py_XDECREF(rate->var_poisson);
@@ -1294,6 +1334,9 @@ create_rateint_product(
  return 0;
 
 FAILED_ALLOC:
+ PyErr_SetString(PyExc_MemoryError, msg);
+ err_ols_print("%s\n", msg);
+
  Py_XDECREF(rateint->slope);
  Py_XDECREF(rateint->dq);
  Py_XDECREF(rateint->var_poisson);
@@ -1606,24 +1649,24 @@ get_ramp_data(
  PyObject * args) /* The C extension module arguments */
 {
  struct ramp_data * rd = calloc(1, sizeof(*rd)); /* Allocate memory */
- PyObject * Py_ramp_data;
+ PyObject * Py_ramp_data = NULL;
  char * msg = "Couldn't allocate memory for ramp data structure.";
 
  /* Make sure memory allocation worked */
  if (NULL==rd) {
  PyErr_SetString(PyExc_MemoryError, msg);
  err_ols_print("%s\n", msg);
- return NULL;
+ goto END;
  }
 
  if (get_ramp_data_parse(&Py_ramp_data, rd, args)) {
  FREE_RAMP_DATA(rd);
- return NULL;
+ goto END;
  }
 
  if (get_ramp_data_arrays(Py_ramp_data, rd)) {
  FREE_RAMP_DATA(rd);
- return NULL;
+ goto END;
  }
 
  /* Void function */
@@ -1643,10 +1686,11 @@ get_ramp_data(
  PyErr_SetString(PyExc_MemoryError, msg);
  err_ols_print("%s\n", msg);
  FREE_RAMP_DATA(rd);
- return NULL;
+ goto END;
  }
  }
 
+END:
  return rd;
 }
 
@@ -1743,6 +1787,8 @@ get_ramp_data_parse(
  return 1;
  }
 
+ /* Note: Freeing 'weight' causes seg fault: 'pointer being freed was not allocated' */
+
  return 0;
 }
 
@@ -2277,6 +2323,7 @@ py_ramp_data_get_float(
 
  Obj = PyObject_GetAttrString(rd, attr);
  val = (float)PyFloat_AsDouble(Obj);
+ Py_XDECREF(Obj);
 
  return val;
 }
@@ -2295,6 +2342,7 @@ py_ramp_data_get_int(
 
  Obj = PyObject_GetAttrString(rd, attr);
  val = (int)PyLong_AsLong(Obj);
+ Py_XDECREF(Obj);
 
  return val;
 }

tests/test_ramp_fitting.py

@@ -1,7 +1,11 @@
 import pytest
 import numpy as np
+
+import sys
+
 from stcal.ramp_fitting.ramp_fit import ramp_fit_data
 from stcal.ramp_fitting.ramp_fit_class import RampData
+from stcal.ramp_fitting.slope_fitter import ols_slope_fitter # c extension
 from stcal.ramp_fitting.utils import compute_num_slices
 
 
@@ -1497,6 +1501,70 @@ def test_one_group():
  assert abs(serr[0, 0] - cerr[0, 0, 0]) < tol
 
 
+def test_compute_num_slices():
+ n_rows = 20
+ max_available_cores = 10
+ assert compute_num_slices("none", n_rows, max_available_cores) == 1
+ assert compute_num_slices("half", n_rows, max_available_cores) == 5
+ assert compute_num_slices("3", n_rows, max_available_cores) == 3
+ assert compute_num_slices("7", n_rows, max_available_cores) == 7
+ assert compute_num_slices("21", n_rows, max_available_cores) == 10
+ assert compute_num_slices("quarter", n_rows, max_available_cores) == 2
+ assert compute_num_slices("7.5", n_rows, max_available_cores) == 1
+ assert compute_num_slices("one", n_rows, max_available_cores) == 1
+ assert compute_num_slices("-5", n_rows, max_available_cores) == 1
+ assert compute_num_slices("all", n_rows, max_available_cores) == 10
+ assert compute_num_slices("3/4", n_rows, max_available_cores) == 1
+ n_rows = 9
+ assert compute_num_slices("21", n_rows, max_available_cores) == 9
+
+
+def test_refcounter():
+ """
+ Get reference of objects before and after C-extension to ensure
+ they are the same.
+ """
+ nints, ngroups, nrows, ncols = 1, 1, 1, 1
+ rnval, gval = 10.0, 5.0
+ frame_time, nframes, groupgap = 10.736, 4, 1
+ # frame_time, nframes, groupgap = 10.736, 1, 0
+
+ dims = nints, ngroups, nrows, ncols
+ var = rnval, gval
+ tm = frame_time, nframes, groupgap
+
+ ramp, gain, rnoise = create_blank_ramp_data(dims, var, tm)
+
+ ramp.data[0, 0, 0, 0] = 105.31459
+
+ b_data = sys.getrefcount(ramp.data)
+ b_dq = sys.getrefcount(ramp.groupdq)
+ b_err = sys.getrefcount(ramp.err)
+ b_pdq = sys.getrefcount(ramp.pixeldq)
+ b_dc = sys.getrefcount(ramp.average_dark_current)
+
+ wt, opt = "optimal", False
+ image, integ, opt= ols_slope_fitter(ramp, gain, rnoise, wt, opt)
+
+ a_data = sys.getrefcount(ramp.data)
+ a_dq = sys.getrefcount(ramp.groupdq)
+ a_err = sys.getrefcount(ramp.err)
+ a_pdq = sys.getrefcount(ramp.pixeldq)
+ a_dc = sys.getrefcount(ramp.average_dark_current)
+
+ # Verify reference counts are not affected by the C-extension, indicating
+ # memory will be properly managed.
+ assert b_data == a_data
+ assert b_dq == a_dq
+ assert b_err == a_err
+ assert b_pdq == a_pdq
+ assert b_dc == a_dc
+
+
+# -----------------------------------------------------------------------------
+# Set up functions
+
+
 def create_blank_ramp_data(dims, var, tm):
  """
  Create empty RampData classes, as well as gain and read noise arrays,
@@ -1531,28 +1599,6 @@ def create_blank_ramp_data(dims, var, tm):
  return ramp_data, gain, rnoise
 
 
-def test_compute_num_slices():
- n_rows = 20
- max_available_cores = 10
- assert compute_num_slices("none", n_rows, max_available_cores) == 1
- assert compute_num_slices("half", n_rows, max_available_cores) == 5
- assert compute_num_slices("3", n_rows, max_available_cores) == 3
- assert compute_num_slices("7", n_rows, max_available_cores) == 7
- assert compute_num_slices("21", n_rows, max_available_cores) == 10
- assert compute_num_slices("quarter", n_rows, max_available_cores) == 2
- assert compute_num_slices("7.5", n_rows, max_available_cores) == 1
- assert compute_num_slices("one", n_rows, max_available_cores) == 1
- assert compute_num_slices("-5", n_rows, max_available_cores) == 1
- assert compute_num_slices("all", n_rows, max_available_cores) == 10
- assert compute_num_slices("3/4", n_rows, max_available_cores) == 1
- n_rows = 9
- assert compute_num_slices("21", n_rows, max_available_cores) == 9
-
-
-# -----------------------------------------------------------------------------
-# Set up functions
-
-
 def setup_inputs(dims, var, tm):
  """
  Given dimensions, variances, and timing data, this creates test data to