Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Cosmics kgain #211

Open
wants to merge 7 commits into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
10 changes: 7 additions & 3 deletions corgidrp/l1_to_l2a.py
Original file line number Diff line number Diff line change
Expand Up @@ -154,7 +154,7 @@ def prescan_biassub(input_dataset, noise_maps=None, return_full_frame=False,

return output_dataset

def detect_cosmic_rays(input_dataset, detector_params, sat_thresh=0.7,
def detect_cosmic_rays(input_dataset, detector_params, k_gain = None, sat_thresh=0.7,
plat_thresh=0.7, cosm_filter=1, cosm_box=3, cosm_tail=10,
mode='image'):
"""
Expand All @@ -166,6 +166,7 @@ def detect_cosmic_rays(input_dataset, detector_params, sat_thresh=0.7,
Args:
input_dataset (corgidrp.data.Dataset): a dataset of Images that need cosmic ray identification (L1-level)
detector_params (corgidrp.data.DetectorParams): a calibration file storing detector calibration values
k_gain (corgidrp.data.KGain): KGain calibration file
sat_thresh (float):
Multiplication factor for the pixel full-well capacity (fwc) that determines saturated cosmic
pixels. Interval 0 to 1, defaults to 0.7. Lower numbers are more aggressive in flagging saturation.
Expand Down Expand Up @@ -205,9 +206,12 @@ def detect_cosmic_rays(input_dataset, detector_params, sat_thresh=0.7,

crmasked_cube = crmasked_dataset.all_data


# Calculate the full well capacity for every frame in the dataset
kgain = np.array([detector_params.params['kgain'] for frame in crmasked_dataset])
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If kgain is None, use the value from detector_params. And k_gain would have to be listed before detector_params in the list of arguments so that the "AUTOMATIC, OPTIONAL" works.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

By the way: is there a reason why an array of equal kgain values is made? We can divide an array by a scalar.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

agree, not really needed to have an array since we are not assuming to deal with different kgain values right now

if k_gain == None:
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

should use if k_gain is None (forget the exact reason, but is is preferred to ==)

kgain = np.array([detector_params.params['kgain'] for frame in crmasked_dataset])
else:
#get the kgain value from the k_gain calibration file
kgain = k_gain.value
emgain_list = []
for frame in crmasked_dataset:
try: # use measured gain if available TODO change hdr name if necessary
Expand Down
3 changes: 2 additions & 1 deletion corgidrp/recipe_templates/build_trad_dark_full_frame.json
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,8 @@
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC"
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

for this and other recipes below, can the indentation be at the same level?

},
"keywords" : {
"mode" : "full"
Expand Down
3 changes: 2 additions & 1 deletion corgidrp/recipe_templates/build_trad_dark_image.json
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,8 @@
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC"
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"
},
"keywords" : {
"mode" : "full"
Expand Down
3 changes: 2 additions & 1 deletion corgidrp/recipe_templates/l1_flat_and_bp.json
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,8 @@
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC"
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"
}
},
{
Expand Down
3 changes: 2 additions & 1 deletion corgidrp/recipe_templates/l1_to_boresight.json
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,8 @@
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC"
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"
}
},
{
Expand Down
7 changes: 7 additions & 0 deletions corgidrp/recipe_templates/l1_to_kgain.json
Original file line number Diff line number Diff line change
Expand Up @@ -31,6 +31,13 @@
}
}
},
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"
}
},
{
"name" : "calibrate_kgain",
"keywords" : {
Expand Down
3 changes: 2 additions & 1 deletion corgidrp/recipe_templates/l1_to_l2a_basic.json
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,8 @@
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC"
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"
}
},
{
Expand Down
3 changes: 2 additions & 1 deletion corgidrp/recipe_templates/l1_to_l2a_eng.json
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,8 @@
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC"
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"
}
},
{
Expand Down
4 changes: 3 additions & 1 deletion corgidrp/recipe_templates/l1_to_l2a_noisemap.json
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,9 @@
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC"
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"

}
},
{
Expand Down
4 changes: 3 additions & 1 deletion corgidrp/recipe_templates/l1_to_l2a_nonlin.json
Original file line number Diff line number Diff line change
Expand Up @@ -34,7 +34,9 @@
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC"
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"

}
},
{
Expand Down
3 changes: 2 additions & 1 deletion corgidrp/recipe_templates/l1_to_l2b.json
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,8 @@
{
"name" : "detect_cosmic_rays",
"calibs" : {
"DetectorParams" : "AUTOMATIC"
"DetectorParams" : "AUTOMATIC",
"KGain" : "AUTOMATIC, OPTIONAL"
}
},
{
Expand Down
9 changes: 8 additions & 1 deletion tests/e2e_tests/l1_to_l2a_e2e.py
Original file line number Diff line number Diff line change
Expand Up @@ -79,14 +79,21 @@ def test_l1_to_l2a(tvacdata_path, e2eoutput_path):
noise_map.save(filedir=l2a_outputdir, filename="mock_detnoisemaps.fits")
this_caldb.create_entry(noise_map)

# KGain
kgain_val = 8.7
kgain = data.KGain(np.array([[kgain_val]]), pri_hdr=pri_hdr, ext_hdr=ext_hdr,
input_dataset=mock_input_dataset)
kgain.save(filedir=l2a_outputdir, filename="mock_kgain.fits")
this_caldb.create_entry(kgain)

####### Run the walker on some test_data

walker.walk_corgidrp(l1_data_filelist, "", l2a_outputdir, template="l1_to_l2a_basic.json")

# clean up by removing entry
this_caldb.remove_entry(nonlinear_cal)
this_caldb.remove_entry(noise_map)

this_caldb.remove_entry(kgain)

##### Check against TVAC data
new_l2a_filenames = [os.path.join(l2a_outputdir, "{0}.fits".format(i)) for i in [90499, 90500]]
Expand Down
12 changes: 10 additions & 2 deletions tests/e2e_tests/nonlin_e2e.py
Original file line number Diff line number Diff line change
Expand Up @@ -152,6 +152,15 @@ def test_nonlin_cal_e2e(
ext_hdr=ext_hdr,
input_dataset=mock_input_dataset)
nonlinear_cal.save(filedir=e2eoutput_path, filename="nonlin_tvac.fits" )


# KGain
kgain_val = 8.7
kgain = data.KGain(np.array([[kgain_val]]), pri_hdr=pri_hdr, ext_hdr=ext_hdr,
input_dataset=mock_input_dataset)
kgain.save(filedir=e2eoutput_path, filename="mock_kgain.fits")
this_caldb = caldb.CalDB()
this_caldb.create_entry(kgain)

# Run the walker on some test_data
print('Running walker')
Expand Down Expand Up @@ -205,9 +214,8 @@ def test_nonlin_cal_e2e(

# remove entry from caldb
nonlin_entry = data.NonLinearityCalibration(os.path.join(e2eoutput_path, nonlin_out_filename))
this_caldb = caldb.CalDB()
this_caldb.remove_entry(nonlin_entry)

this_caldb.remove_entry(kgain)
# Print success message
print('e2e test for NL passed')

Expand Down
26 changes: 15 additions & 11 deletions tests/test_cr_detection.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,11 @@
non_linearity_correction = data.NonLinearityCalibration(tvac_nonlin_data,pri_hdr=pri_hdr,ext_hdr=ext_hdr,input_dataset = dummy_dataset)
non_linearity_correction.save(filename = nonlin_fits_filepath)

# Make a dummy kgain calibration file
kgain = 8.7
gain_value = np.array([[kgain]])
k_gain = data.KGain(gain_value, pri_hdr = pri_hdr, ext_hdr = ext_hdr, input_dataset = dummy_dataset)

## Copy-pasted II&T code from https://github.com/roman-corgi/cgi_iit_drp/blob/main/proc_cgi_frame_NTR/proc_cgi_frame/gsw_remove_cosmics.py ##

def find_plateaus_iit(streak_row, fwc, sat_thresh, plat_thresh, cosm_filter):
Expand Down Expand Up @@ -198,7 +203,6 @@ def test_iit_vs_corgidrp():
and check that output is consistent with results II&T code.
"""

kgain = detector_params.params['kgain']
fwc_em = detector_params.params['fwc_em'] / kgain
fwc_pp = detector_params.params['fwc_pp'] / kgain
em_gain = 500
Expand Down Expand Up @@ -233,7 +237,7 @@ def test_iit_vs_corgidrp():
iit_masks_arr = np.array(iit_masks)

# corgidrp version
crmasked_dataset = detect_cosmic_rays(dataset, detector_params, sat_thresh,
crmasked_dataset = detect_cosmic_rays(dataset, detector_params, k_gain, sat_thresh,
plat_thresh, cosm_filter, cosm_box,
cosm_tail, mode)
corgi_crmask_bool = np.where(crmasked_dataset.all_dq>0,1,0)
Expand Down Expand Up @@ -264,7 +268,7 @@ def test_correct_headers():
"""
# create simulated data
dataset = mocks.create_cr_dataset(nonlin_fits_filepath, filedir=datadir, numfiles=2,numCRs=5, plateau_length=10)
output_dataset = detect_cosmic_rays(dataset, detector_params)
output_dataset = detect_cosmic_rays(dataset, detector_params, k_gain)

for frame in output_dataset:
if not ("FWC_EM_E" in frame.ext_hdr):
Expand Down Expand Up @@ -377,7 +381,7 @@ def test_mask():
check_mask = np.zeros_like(dataset.all_dq, dtype=int)
c_tail = 6
check_mask[0,1, 2:2+cosm_filter+c_tail+1] = 1 #add 1 to include last column in the slice
dataset_masked = detect_cosmic_rays(dataset, detector_params, sat_thresh, plat_thresh, cosm_filter, cosm_box=0, cosm_tail=c_tail)
dataset_masked = detect_cosmic_rays(dataset, detector_params, k_gain, sat_thresh, plat_thresh, cosm_filter, cosm_box=0, cosm_tail=c_tail)
if not np.where(dataset_masked.all_dq>0,1,0) == approx(check_mask):
raise Exception("Incorrect pixels were masked.")

Expand All @@ -402,14 +406,14 @@ def test_mask_box():
frame = data.Image(bs_image_box, pri_hdr=prihdr,
ext_hdr=exthdr)
dataset = data.Dataset([frame])
dataset_masked = detect_cosmic_rays(dataset, detector_params, sat_thresh,
dataset_masked = detect_cosmic_rays(dataset, detector_params, k_gain, sat_thresh,
plat_thresh, cosm_filter=2, cosm_box=0,
cosm_tail=20)

assert not (np.array_equal(np.where(dataset_masked.all_dq>0,1,0)[0], check_mask)) # since cosm_box=0

# now use cosm_box=2 to catch pixels surrounding head
dataset_masked2 = detect_cosmic_rays(dataset, detector_params, sat_thresh,
dataset_masked2 = detect_cosmic_rays(dataset, detector_params, k_gain, sat_thresh,
plat_thresh, cosm_filter=2, cosm_box=2,
cosm_tail=20)

Expand Down Expand Up @@ -442,7 +446,7 @@ def test_mask_box_corners():
frame = data.Image(image, pri_hdr=prihdr,
ext_hdr=exthdr)
dataset = data.Dataset([frame])
dataset_masked = detect_cosmic_rays(dataset, detector_params, sat_thresh,
dataset_masked = detect_cosmic_rays(dataset, detector_params, k_gain, sat_thresh,
plat_thresh, cosm_filter=2, cosm_box=2)

if not np.array_equal(np.where(dataset_masked.all_dq>0,1,0)[0], check_mask):
Expand All @@ -467,7 +471,7 @@ def test_cosm_tail_2():
frame = data.Image(image, pri_hdr=prihdr,
ext_hdr=exthdr)
dataset = data.Dataset([frame])
dataset_masked = detect_cosmic_rays(dataset, detector_params, sat_thresh,
dataset_masked = detect_cosmic_rays(dataset, detector_params, k_gain, sat_thresh,
plat_thresh, cosm_filter=2, cosm_box=0,
cosm_tail=1)

Expand All @@ -479,7 +483,7 @@ def test_cosm_tail_2():
# cosmic head #2
check_mask[-2,6:6+2+3+1] = 1

dataset_masked = detect_cosmic_rays(dataset, detector_params, sat_thresh,
dataset_masked = detect_cosmic_rays(dataset, detector_params, k_gain, sat_thresh,
plat_thresh, cosm_filter=2, cosm_box=0,
cosm_tail=3)

Expand All @@ -501,7 +505,7 @@ def test_cosm_tail_bleed_over():
frame = data.Image(image, pri_hdr=prihdr,
ext_hdr=exthdr)
dataset = data.Dataset([frame])
dataset_masked = detect_cosmic_rays(dataset, detector_params, sat_thresh,
dataset_masked = detect_cosmic_rays(dataset, detector_params, k_gain, sat_thresh,
plat_thresh, cosm_filter=2, cosm_box=2,
cosm_tail=14, mode='full')

Expand All @@ -511,7 +515,7 @@ def test_cosm_tail_bleed_over():
check_mask[-1,0:12] = 0 # undo the bleed over
# cosm_box=2 again since I undid some in previous line
check_mask[-4:,4:9] = 1
dataset_masked = detect_cosmic_rays(dataset, detector_params, sat_thresh,
dataset_masked = detect_cosmic_rays(dataset, detector_params, k_gain, sat_thresh,
plat_thresh, cosm_filter=2, cosm_box=2,
cosm_tail=14)

Expand Down
21 changes: 21 additions & 0 deletions tests/test_walker.py
Original file line number Diff line number Diff line change
Expand Up @@ -63,6 +63,16 @@ def test_autoreducing():
new_nonlinearity.ext_hdr.set('DRPVERSN', corgidrp.__version__, "corgidrp version that produced this file")
mycaldb = caldb.CalDB()
mycaldb.create_entry(new_nonlinearity)

#Make a KGain calibration file
kgain_arr = np.array([[8.8]])
new_kgain = data.KGain(kgain_arr,pri_hdr=prihdr,ext_hdr=exthdr,input_dataset = dummy_dataset)
new_kgain.ext_hdr.set('DRPCTIME', time.Time.now().isot, "When this file was saved")
new_kgain.ext_hdr.set('DRPVERSN', corgidrp.__version__, "corgidrp version that produced this file")
new_kgain.save(filedir = os.path.join(os.path.dirname(__file__), "test_data"), filename = "kgain.fits")

mycaldb.create_entry(new_kgain)


CPGS_XML_filepath = "" # not yet implemented

Expand All @@ -85,6 +95,7 @@ def test_autoreducing():

# clean up
mycaldb.remove_entry(new_nonlinearity)
mycaldb.remove_entry(new_kgain)

def test_auto_template_identification():
"""
Expand Down Expand Up @@ -432,6 +443,15 @@ def test_jit_calibs():
mycaldb = caldb.CalDB()
mycaldb.create_entry(new_nonlinearity)

#Make a KGain calibration file
kgain_arr = np.array([[8.8]])
new_kgain = data.KGain(kgain_arr,pri_hdr=prihdr,ext_hdr=exthdr,input_dataset = dummy_dataset)
new_kgain.ext_hdr.set('DRPCTIME', time.Time.now().isot, "When this file was saved")
new_kgain.ext_hdr.set('DRPVERSN', corgidrp.__version__, "corgidrp version that produced this file")
new_kgain.save(filedir = os.path.join(os.path.dirname(__file__), "test_data"), filename = "kgain.fits")

mycaldb.create_entry(new_kgain)

CPGS_XML_filepath = "" # not yet implemented

# generate recipe and check we haven't defined anyhting yet
Expand Down Expand Up @@ -494,6 +514,7 @@ def test_jit_calibs():

# clean up
mycaldb.remove_entry(new_nonlinearity)
mycaldb.remove_entry(new_kgain)

corgidrp.jit_calib_id = old_setting

Expand Down
Loading