diff --git a/CHANGES.rst b/CHANGES.rst index 31bed2da..2f3a46a9 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -16,6 +16,11 @@ Bug Fixes - +ramp_fitting +~~~~~~~~~~~~ + +- Fixed memory leak in C-extension.[#281] + 1.8.0 (2024-08-14) ================== diff --git a/src/stcal/ramp_fitting/src/slope_fitter.c b/src/stcal/ramp_fitting/src/slope_fitter.c index b27a2f35..860c601e 100644 --- a/src/stcal/ramp_fitting/src/slope_fitter.c +++ b/src/stcal/ramp_fitting/src/slope_fitter.c @@ -3,9 +3,13 @@ #include #include -#include #include +#include +#include #include +#include +#include +#include #include #include @@ -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; } diff --git a/tests/test_ramp_fitting.py b/tests/test_ramp_fitting.py index 3cdd2215..412c6ab8 100644 --- a/tests/test_ramp_fitting.py +++ b/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