Skip to content

Commit

Permalink
Updating debugging functions for the RampData class and updating the …
Browse files Browse the repository at this point in the history 
…C median rate computation, fixing many differences with the python code for MIRI image regression test.
  • Loading branch information
@kmacdonald-stsci
kmacdonald-stsci committed Feb 2, 2024
1 parent f3d07e4 commit 855b27b
Show file tree
Hide file tree
Showing 2 changed files with 95 additions and 58 deletions.
133 changes: 80 additions & 53 deletions src/stcal/ramp_fitting/ramp_fit_class.py
View file
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
INDENT = " "

class RampData:
def __init__(self):
"""Creates an internal ramp fit class."""
Expand Down Expand Up @@ -191,73 +193,98 @@ def dbg_print_pixel_info(self, row, col):
# print(f" err :\n{self.err[:, :, row, col]}")
# print(f" pixeldq :\n{self.pixeldq[row, col]}")

def dbg_write_ramp_data_pix(self, fname, row, col, gain, rnoise):
print(f"*** {fname} ***")
with open(fname, "w") as fd:
fd.write("def create_ramp_data_pixel():\n")
indent = " "
fd.write(f"{indent}'''\n")
fd.write(f"{indent}Using pixel ({row}, {col}).n")
fd.write(f"{indent}'''\n")
fd.write(f"{indent}ramp_data = RampData()\n\n")
def dbg_write_ramp_data_pix_pre(self, fname, row, col, fd):
fd.write("def create_ramp_data_pixel():\n")
indent = INDENT
fd.write(f"{indent}'''\n")
fd.write(f"{indent}Using pixel ({row}, {col}).n")
fd.write(f"{indent}'''\n")
fd.write(f"{indent}ramp_data = RampData()\n\n")

fd.write(f"{indent}ramp_data.instrument_name = {self.instrument_name}\n\n")

fd.write(f"{indent}ramp_data.frame_time = {self.frame_time}\n")
fd.write(f"{indent}ramp_data.group_time = {self.group_time}\n")
fd.write(f"{indent}ramp_data.groupgap = {self.groupgap}\n")
fd.write(f"{indent}ramp_data.nframes = {self.nframes}\n")
fd.write(f"{indent}ramp_data.drop_frames1 = {self.drop_frames1}\n\n")

fd.write(f"{indent}ramp_data.flags_do_not_use = {self.flags_do_not_use}\n")
fd.write(f"{indent}ramp_data.flags_jump_det = {self.flags_jump_det}\n")
fd.write(f"{indent}ramp_data.flags_saturated = {self.flags_saturated}\n")
fd.write(f"{indent}ramp_data.flags_no_gain_val = {self.flags_no_gain_val}\n")
fd.write(f"{indent}ramp_data.flags_unreliable_slope = {self.flags_unreliable_slope}\n\n")

fd.write(f"{indent}ramp_data.instrument_name = {self.instrument_name}\n\n")

fd.write(f"{indent}ramp_data.frame_time = {self.frame_time}\n")
fd.write(f"{indent}ramp_data.group_time = {self.group_time}\n")
fd.write(f"{indent}ramp_data.groupgap = {self.groupgap}\n")
fd.write(f"{indent}ramp_data.nframes = {self.nframes}\n")
fd.write(f"{indent}ramp_data.drop_frames1 = {self.drop_frames1}\n\n")
fd.write(f"{indent}ramp_data.start_row = 0\n")
fd.write(f"{indent}ramp_data.num_rows = 1\n\n")

fd.write(f"{indent}ramp_data.flags_do_not_use = {self.flags_do_not_use}\n")
fd.write(f"{indent}ramp_data.flags_jump_det = {self.flags_jump_det}\n")
fd.write(f"{indent}ramp_data.flags_saturated = {self.flags_saturated}\n")
fd.write(f"{indent}ramp_data.flags_no_gain_val = {self.flags_no_gain_val}\n")
fd.write(f"{indent}ramp_data.flags_unreliable_slope = {self.flags_unreliable_slope}\n\n")
fd.write(f"{indent}ramp_data.suppress_one_group_ramps = {self.suppress_one_group_ramps}\n\n")

nints, ngroups, nrows, ncols = self.data.shape
fd.write(f"{indent}data = np.zeros(({nints}, {ngroups}, 1, 1), dtype=np.float32)\n")
fd.write(f"{indent}err = np.zeros(({nints}, {ngroups}, 1, 1), dtype=np.float32)\n")
fd.write(f"{indent}gdq = np.zeros(({nints}, {ngroups}, 1, 1), dtype=np.uint8)\n")
fd.write(f"{indent}pdq = np.zeros((1, 1), dtype=np.uint32)\n")

fd.write(f"{indent}ramp_data.start_row = 0\n")
fd.write(f"{indent}ramp_data.num_rows = 1\n\n")

fd.write(f"{indent}ramp_data.suppress_one_group_ramps = {self.suppress_one_group_ramps}\n\n")
def dbg_write_ramp_data_pix_post(self, fname, row, col, fd):
indent = INDENT

nints, ngroups, nrows, ncols = self.data.shape
fd.write(f"{indent}data = np.zeros(({nints}, {ngroups}, 1, 1), dtype=np.float32)\n")
fd.write(f"{indent}err = np.zeros(({nints}, {ngroups}, 1, 1), dtype=np.float32)\n")
fd.write(f"{indent}gdq = np.zeros(({nints}, {ngroups}, 1, 1), dtype=np.uint8)\n")
fd.write(f"{indent}pdq = np.zeros((1, 1), dtype=np.uint32)\n")
fd.write(f"{indent}ramp_data.data = data\n")
fd.write(f"{indent}ramp_data.err = err\n")
fd.write(f"{indent}ramp_data.groupdq = gdq\n")
fd.write(f"{indent}ramp_data.pixeldq = pdq\n")
fd.write(f"{indent}ramp_data.zeroframe = zframe\n\n")

import numpy as np
fd.write(f"{indent}return ramp_data, ngain, nrnoise\n")

arr_str = np.array2string(self.data[:, :, row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}data[:, :, 0, 0] = np.array({arr_str})\n\n")
arr_str = np.array2string(self.err[:, :, row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}err[:, :, 0, 0] = np.array({arr_str})\n\n")
arr_str = np.array2string(self.groupdq[:, :, row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}gdq[:, :, 0, 0] = np.array({arr_str})\n\n")
arr_str = np.array2string(self.pixeldq[row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}pdq[0, 0] = np.array({arr_str})\n\n")
def dbg_write_ramp_data_pix_pixel(self, fname, row, col, gain, rnoise, fd):
import numpy as np
indent = INDENT

if self.zeroframe is not None:
fd.write(f"{indent}zframe = np.zeros((1, 1), dtype=np.float32)\n\n")
arr_str = np.array2string(self.zeroframe[row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}zframe[0, 0] = np.array({arr_str})\n\n")
else:
fd.write(f"{indent}zframe = None\n\n")
# XXX Make this a separate function
delimiter = "-" * 40
fd.write(f"{indent}# {delimiter}\n\n");
fd.write(f"{indent}# ({row}, {col})\n\n");

fd.write(f"{indent}ramp_data.data = data\n")
fd.write(f"{indent}ramp_data.err = err\n")
fd.write(f"{indent}ramp_data.groupdq = gdq\n")
fd.write(f"{indent}ramp_data.pixeldq = pdq\n")
fd.write(f"{indent}ramp_data.zeroframe = zframe\n\n")
nints = self.data.shape[0]

fd.write(f"{indent}ngain = np.zeros((1, 1), dtype=np.float32)\n\n")
fd.write(f"{indent}ngain[0, 0] = np.array([[{gain[row, col]}]])\n\n")
for integ in range(nints):
arr_str = np.array2string(self.data[integ, :, row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}data[{integ}, :, 0, 0] = np.array({arr_str})\n")
fd.write("\n")

fd.write(f"{indent}nrnoise = np.zeros((1, 1), dtype=np.float32)\n\n")
fd.write(f"{indent}nrnoise[0, 0] = np.array([[{rnoise[row, col]}]])\n\n")
for integ in range(nints):
arr_str = np.array2string(self.err[integ, :, row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}err[{integ}, :, 0, 0] = np.array({arr_str})\n")
fd.write("\n")

fd.write(f"{indent}return ramp_data, ngain, nrnoise\n")
for integ in range(nints):
arr_str = np.array2string(self.groupdq[integ, :, row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}gdq[{integ}, :, 0, 0] = np.array({arr_str})\n")
fd.write("\n")

arr_str = np.array2string(self.pixeldq[row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}pdq[0, 0] = {arr_str}\n\n")

if self.zeroframe is not None:
fd.write(f"{indent}zframe = np.zeros((1, 1), dtype=np.float32)\n\n")
arr_str = np.array2string(self.zeroframe[row, col], precision=12, max_line_width=np.nan, separator=", ")
fd.write(f"{indent}zframe[0, 0] = {arr_str}\n\n")
else:
fd.write(f"{indent}zframe = None\n\n")

fd.write(f"{indent}ngain = np.zeros((1, 1), dtype=np.float32)\n")
fd.write(f"{indent}ngain[0, 0] = {gain[row, col]}\n\n")

fd.write(f"{indent}nrnoise = np.zeros((1, 1), dtype=np.float32)\n")
fd.write(f"{indent}nrnoise[0, 0] = {rnoise[row, col]}\n\n")


def dbg_write_ramp_data_pix(self, fname, row, col, gain, rnoise):
print(f"*** {fname} ***")
with open(fname, "w") as fd:
self.dbg_write_ramp_data_pix_pre(fname, row, col, fd)
self.dbg_write_ramp_data_pix_pixel(fname, row, col, gain, rnoise, fd)
self.dbg_write_ramp_data_pix_post(fname, row, col, fd)
20 changes: 15 additions & 5 deletions src/stcal/ramp_fitting/src/slope_fitter.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -461,7 +461,7 @@ median_rate_integration(
real_t * mrate, real_t * int_data, uint8_t * int_dq,
struct ramp_data * rd, struct pixel_ramp * pr);

static void
static int
median_rate_integration_sort(
real_t * loc_integ, uint8_t * int_dq,
struct ramp_data * rd, struct pixel_ramp * pr);
Expand Down Expand Up @@ -1973,6 +1973,7 @@ median_rate_integration(
real_t * loc_integ = (real_t*)calloc(pr->ngroups, sizeof(*loc_integ));
const char * msg = "Couldn't allocate memory for integration median rate.";
npy_intp k, loc_integ_len;
int nan_cnt;

/* Make sure memory allocation worked */
if (NULL==loc_integ) {
Expand All @@ -1992,7 +1993,7 @@ median_rate_integration(
}

/* Sort first differences with NaN's based on DQ flags */
median_rate_integration_sort(loc_integ, int_dq, rd, pr);
nan_cnt = median_rate_integration_sort(loc_integ, int_dq, rd, pr);

/*
* Get the NaN median using the sorted first differences. Note that the
Expand All @@ -2015,7 +2016,7 @@ median_rate_integration(
* Set flagged groups to NaN.
* Sort the modified data.
*/
static void
static int
median_rate_integration_sort(
real_t * loc_integ,
uint8_t * int_dq,
Expand All @@ -2024,10 +2025,11 @@ median_rate_integration_sort(
{
npy_intp k, ngroups = pr->ngroups;
real_t loc0 = loc_integ[0];
int nan_cnt = 0;

/* Compute first differences */
if (1 == ngroups ) {
return;
return nan_cnt;
} else {
for (k=0; k<ngroups-1; ++k) {
if (rd->jump & int_dq[k+1]) {
Expand All @@ -2036,14 +2038,21 @@ median_rate_integration_sort(
} else {
loc_integ[k] = loc_integ[k+1] - loc_integ[k];
}
if (isnan(loc_integ[k])) {
nan_cnt++;
}
}
}
if (isnan(loc_integ[0])) {

// if (isnan(loc_integ[0])) {
if (nan_cnt == (ngroups-1)) {
loc_integ[0] = loc0;
}

/* NaN sort first differences */
qsort(loc_integ, ngroups-1, sizeof(loc_integ[0]), median_rate_integration_sort_cmp);

return nan_cnt;
}

/* The comparison function for qsort with NaN's */
Expand Down Expand Up @@ -2270,6 +2279,7 @@ ramp_fit_pixel(
print_delim();
dbg_ols_print("Pixel (%ld, %ld)\n", pr->row, pr->col);
dbg_ols_print("Median Rate = %.10f\n", pr->median_rate);
print_delim();
#endif

/* Clean up any thing from the last pixel ramp */
Expand Down

0 comments on commit 855b27b

Please sign in to comment.