From f92ac6ea740bf83028b24fe7dd31d4ee7ddeb2a6 Mon Sep 17 00:00:00 2001 From: Philip Hodge Date: Thu, 31 Oct 2019 13:13:43 -0400 Subject: [PATCH 1/6] Add unit tests for defringe Unit tests were added for defringe.py. Changes were made to defringe.py to print more information with regard to zero or negative values in the fringe flat and data quality flags (if any) in the fringe flat. --- stistools/defringe.py | 20 ++-- tests/test_defringe.py | 217 +++++++++++++++++++++++++++++++++++++++++ 2 files changed, 231 insertions(+), 6 deletions(-) create mode 100644 tests/test_defringe.py diff --git a/stistools/defringe.py b/stistools/defringe.py index 83dc2d9d..c4ae0eff 100644 --- a/stistools/defringe.py +++ b/stistools/defringe.py @@ -85,18 +85,25 @@ def defringe(science_file, fringe_flat, overwrite=True, verbose=True): # Since we're going to divide by fringe_data, make sure there aren't any # pixels where it's zero. There shouldn't be any negative values, either. fringe_mask = (fringe_data <= 0.) + n_not_positive = fringe_mask.sum() + if n_not_positive > 0 and verbose: + print('{} pixels in the fringe flat were less than or equal to 0' + .format(n_not_positive)) if fringe_dq is not None: temp = np.bitwise_and(fringe_dq, sdqflags) fringe_dq_mask = (temp > 0) del temp - fringe_mask = np.logical_or(fringe_mask, fringe_dq_mask) + n_fringe_dq = fringe_dq_mask.sum() + if n_fringe_dq > 0: + fringe_mask = np.logical_or(fringe_mask, fringe_dq_mask) + if verbose: + print('{} pixels were flagged as bad in the fringe flat DQ array' + .format(n_fringe_dq)) # For pixels that are bad in the fringe flat, we will not make any change # to the science data. + # Note: Save fringe_mask and n_fringe_mask for later. n_fringe_mask = fringe_mask.sum() if n_fringe_mask > 0: - if verbose: - print('{} pixels in the fringe flat were less than or equal to 0' - .format(n_fringe_mask)) fringe_data[fringe_mask] = 1. # Correct the data in the science file: @@ -140,7 +147,8 @@ def defringe(science_file, fringe_flat, overwrite=True, verbose=True): if science_dq is None: science_dq = np.zeros(science_data.shape, dtype=np.int16) if n_fringe_mask > 0: - # Flag pixels that had fringe flat data <= 0. + # Flag pixels that had fringe flat data <= 0 or were flagged + # as bad in the fringe flat DQ array. science_dq[fringe_mask] |= 512 # bad pixel in ref file if fringe_dq is not None: # Combine fringe flat DQ with science DQ. @@ -163,7 +171,7 @@ def defringe(science_file, fringe_flat, overwrite=True, verbose=True): print('Removing and recreating {}'.format(drj_filename)) os.remove(drj_filename) science_hdu.writeto(drj_filename) - print('Defringed science saved to {}'.format(drj_filename)) + print('Defringed science data saved to {}'.format(drj_filename)) return drj_filename diff --git a/tests/test_defringe.py b/tests/test_defringe.py new file mode 100644 index 00000000..36487e9d --- /dev/null +++ b/tests/test_defringe.py @@ -0,0 +1,217 @@ +import os +import tempfile +import numpy as np + +from astropy.io import fits + +from stistools.defringe import defringe + +shape = (5, 7) +datasize = shape[0] * shape[1] + +# Name of a temporary directory. See functions save_file and clean_up. +temp_dir = None +# A list of the names (including directory) of temporary files created +# for these tests. +filelist = [] # list of names (with directory) of temporary files + + +def save_file_name(filename): + + global filelist, temp_dir + + if temp_dir is None: + temp_dir = tempfile.mkdtemp() + + if len(filename) > len(os.path.basename(filename)): + # `filename` already includes the directory name + fullname = filename + else: + # `filename` does not include the directory yet + fullname = os.path.join(temp_dir, filename) + + # If the current file is already in `filelist`, don't include it again. + if fullname not in filelist: + filelist.append(fullname) + + return fullname + + +def save_file(fd, filename): + + fullname = save_file_name(filename) + + fd.writeto(fullname, overwrite=True) + + return fullname + + +def clean_up(): + + global filelist, temp_dir + + if temp_dir is None: + return + + while filelist: + fullname = filelist.pop() + os.remove(fullname) + os.rmdir(temp_dir) + temp_dir = None + + +def mk_sci_file(filename, n_imsets=None): + + fd = fits.HDUList(fits.PrimaryHDU()) + + dqflags = [4, 8, 512, 0] + if n_imsets is None: + n_imsets = 1 + for extver in range(1, n_imsets + 1): + # SCI HDU + if extver == 3: + hdu = fits.ImageHDU(name="SCI", ver=extver) # no data + else: + sci = np.arange(datasize, dtype=np.float32).reshape(shape) + sci += extver + hdu = fits.ImageHDU(data=sci.copy(), name="SCI", ver=extver) + fd.append(hdu) + # DQ HDU + if extver == 2: + hdu = fits.ImageHDU(name="DQ", ver=extver) # no data + else: + i = extver % 4 + dq = np.zeros(shape, dtype=np.int16) + dq[i, i] = dqflags[i] + hdu = fits.ImageHDU(data=dq.copy(), name="DQ", ver=extver) + fd.append(hdu) + # ERR HDU + if extver == 1: + hdu = fits.ImageHDU(name="ERR", ver=extver) # no data + else: + err = np.ones(shape, dtype=np.float32) + hdu = fits.ImageHDU(data=err.copy(), name="ERR", ver=extver) + fd.append(hdu) + + fullname = save_file(fd, filename) + fd.close() + + return fullname + + +def mk_fringe_flat(filename, hdunum=1): + + if hdunum not in [0, 1]: + print('Warning (mk_fringe_flat): hdunum must be 0 or 1') + + if hdunum == 0: + data = np.ones(shape, dtype=np.float32) * 2. + fd = fits.HDUList(fits.PrimaryHDU(data=data)) + else: + extver = 1 + fd = fits.HDUList(fits.PrimaryHDU()) + data = np.ones(shape, dtype=np.float32) * 2. + data[1, 4] = -1. + data[1, 5] = 0. + hdu = fits.ImageHDU(data=data.copy(), name="SCI", ver=extver) + fd.append(hdu) + dq = np.zeros(shape, dtype=np.int16) + dq[2, 2:4] = 15 + hdu = fits.ImageHDU(data=dq.copy(), name="DQ", ver=extver) + fd.append(hdu) + + fullname = save_file(fd, filename) + fd.close() + + return fullname + + +def test_defringe(): + + n_imsets = 5 + + science_file = mk_sci_file('ttt_crj.fits', n_imsets=n_imsets) + + fringe_flat = mk_fringe_flat('fringe_flat.fits', hdunum=1) + + drj_filename = defringe(science_file, fringe_flat, + overwrite=True, verbose=True) + drj_filename = save_file_name(drj_filename) + + in_sci = fits.open(science_file) + in_flat = fits.open(fringe_flat) + out_sci = fits.open(drj_filename) + + fringe_data = in_flat[('sci', 1)].data + fringe_dq = in_flat[('dq', 1)].data + + results = [] + + # Test 1 + check = len(out_sci) == len(in_sci) + if not check: + print('FAILED: the input and output science files do not have the same number of HDUs') + results.append(check) + + # Test 2 + # Check that the science data array was divided by the fringe flat, + # except at pixels where the fringe flat was less than or equal to zero, + # or the fringe flat was flagged as bad in its DQ array. + f_flat = fringe_data.copy() # because we're going to modify it + sci5_data = in_sci[('sci', 5)].data + ff_is_bad = (fringe_data <= 0.) + f_flat[ff_is_bad] = 1. + dq_flagged = (fringe_dq > 0) + f_flat[dq_flagged] = 1. + compare = sci5_data / f_flat + check = np.allclose(out_sci[('sci', 5)].data, compare, rtol=1.e-6) + if not check: + print('FAILED: science array not flat fielded correctly') + results.append(check) + + # Test 3 + f_flat = fringe_data.copy() + err4_data = in_sci[('err', 4)].data + ff_is_bad = (fringe_data <= 0.) + f_flat[ff_is_bad] = 1. + dq_flagged = (fringe_dq > 0) + f_flat[dq_flagged] = 1. + compare = err4_data / f_flat + check = np.allclose(out_sci[('err', 4)].data, compare, rtol=1.e-6) + if not check: + print('FAILED: error array not flat fielded correctly') + results.append(check) + + # Test 4 + # The DQ extension in imset 2 in the input science file did not have a + # data block, but one should have been created and added to this imset + # in the output science file. + dq2_data = out_sci[('dq', 2)].data + if dq2_data is None: + print("FAILED: out_sci[('dq', 2)].data is None") + results.append(dq2_data is not None) + + # Test 5 + dq4_data = out_sci[('dq', 4)].data + # Pixels [1, 4], [1, 5], [2, 2], and [2, 3] should all be flagged with 512. + # The first two were negative or zero in the fringe flat data, and the + # latter two were flagged with 15 (4 and 8 are both "serious") in the + # fringe flat data quality array. + # The input science DQ array had value 4 at pixel [0, 0] for imset 4, and + # that value should be propagated to the output DQ array. + compare = np.zeros(shape, dtype=np.int16) + compare[1, 4] = 512 + compare[1, 5] = 512 + compare[2, 2] = 512 | 15 + compare[2, 3] = 512 | 15 + compare[0, 0] = 4 + diff = dq4_data - compare + if len(diff.nonzero()[0]) > 0: + print("FAILED: out_sci[('dq', 4)].data is not what was expected") + results.append(len(diff.nonzero()[0]) == 0) + + in_sci.close() + in_flat.close() + in_sci.close() + + clean_up() # remove temporary files and directory From 2854bf9adf9c991fb711d36e4425b2e76e2ebd44 Mon Sep 17 00:00:00 2001 From: Philip Hodge Date: Thu, 31 Oct 2019 15:54:47 -0400 Subject: [PATCH 2/6] Add an assert statement for the tests. --- tests/test_defringe.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/tests/test_defringe.py b/tests/test_defringe.py index 36487e9d..45817f84 100644 --- a/tests/test_defringe.py +++ b/tests/test_defringe.py @@ -210,6 +210,8 @@ def test_defringe(): print("FAILED: out_sci[('dq', 4)].data is not what was expected") results.append(len(diff.nonzero()[0]) == 0) + assert np.alltrue(results) + in_sci.close() in_flat.close() in_sci.close() From 3706c401b3c0e02c8f776f49ef65c516aedc5498 Mon Sep 17 00:00:00 2001 From: Philip Hodge Date: Fri, 1 Nov 2019 10:35:35 -0400 Subject: [PATCH 3/6] Move the assert statement to a point after deleting the temporary files --- tests/test_defringe.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/test_defringe.py b/tests/test_defringe.py index 45817f84..b2ef98b9 100644 --- a/tests/test_defringe.py +++ b/tests/test_defringe.py @@ -210,10 +210,10 @@ def test_defringe(): print("FAILED: out_sci[('dq', 4)].data is not what was expected") results.append(len(diff.nonzero()[0]) == 0) - assert np.alltrue(results) - in_sci.close() in_flat.close() in_sci.close() clean_up() # remove temporary files and directory + + assert np.alltrue(results) From dc3160b1897476994966ed174a4e81f560bc668f Mon Sep 17 00:00:00 2001 From: Philip Hodge Date: Fri, 1 Nov 2019 10:57:22 -0400 Subject: [PATCH 4/6] Remove a line that is no longer needed. --- stistools/defringe.py | 1 - 1 file changed, 1 deletion(-) diff --git a/stistools/defringe.py b/stistools/defringe.py index c4ae0eff..3a12243f 100644 --- a/stistools/defringe.py +++ b/stistools/defringe.py @@ -9,7 +9,6 @@ import numpy as np -stistools.r_util.NOT_APPLICABLE = 'n/a' # Fix a bug in expandFileName() from stistools.r_util import expandFileName # 4 bad detector pixel or beyond aperture From 0d06b412b9425f179a4f0a35b17ea1e03de7e478 Mon Sep 17 00:00:00 2001 From: Philip Hodge Date: Fri, 1 Nov 2019 14:34:24 -0400 Subject: [PATCH 5/6] Add output argument; use "_defringe" as default output suffix --- stistools/defringe.py | 50 ++++++++++++++++++++++++++----------------- 1 file changed, 30 insertions(+), 20 deletions(-) mode change 100644 => 100755 stistools/defringe.py diff --git a/stistools/defringe.py b/stistools/defringe.py old mode 100644 new mode 100755 index 3a12243f..034fe44a --- a/stistools/defringe.py +++ b/stistools/defringe.py @@ -17,35 +17,41 @@ sdqflags = 4 + 8 + 512 -def defringe(science_file, fringe_flat, overwrite=True, verbose=True): +def defringe(science_file, fringe_flat, output=None, + overwrite=True, verbose=True): # dark_file=None, pixel_flat=None """Defringe by dividing the science spectrum by the fringe flat Parameters ---------- science_file : str - The name of the input science file. + The name (including directory) of the input science file. fringe_flat : str - The name of the input fringe flat file. This is the output from - `mkfringeflat`. + The name (including directory) of the input fringe flat file. This + is the output from `mkfringeflat`. + + output : str or None + If `output` is None, the name of the output file will be constructed + by appending "_defringe" to the root of the input file name. + If `output` is not None, this is the name (not including directory) + of the output file. The output will be written to the same directory + as the input science file. overwrite : bool - The name of the output file will be constructed from the name of the - input science file (`science_file`) by replacing the suffix with - 'drj'. If the input name has suffix 'drj', a RuntimeError will be - raised, rather than modifying the input in-place. If there is an existing file with the same name as the output name, the existing file will be overwritten if `overwrite` is True (the default is True). + If the input and output names are the same, a RuntimeError will be + raised, rather than modifying the input in-place. verbose : bool If True (the default), print more info. Returns ------- - drj_filename : str - The name of the output file. + output_filename : str + The name of the output file, including directory. """ if science_file.endswith("_raw.fits"): @@ -55,9 +61,13 @@ def defringe(science_file, fringe_flat, overwrite=True, verbose=True): # Define new filename: science_file = os.path.normpath(expandFileName(science_file)) # Expand IRAF and UNIX $VARS sci_dir, sci_filename = os.path.split(science_file) - sci_root = re.split('\.fits.*', sci_filename, flags=re.IGNORECASE)[0].rsplit('_',1)[0] - drj_filename = os.path.join(sci_dir, sci_root + '_drj.fits') - if science_file == drj_filename: + if output is None: + sci_root = re.split('\.fits.*', sci_filename, + flags=re.IGNORECASE)[0].rsplit('_',1)[0] + output_filename = os.path.join(sci_dir, sci_root + '_defringe.fits') + else: + output_filename = os.path.join(sci_dir, output) + if science_file == output_filename: raise RuntimeError('The input and output file names cannot be the same.') # Get the data from the fringe flat file: @@ -166,13 +176,13 @@ def defringe(science_file, fringe_flat, overwrite=True, verbose=True): # Write to a new FITS file # Remove old version, if it exists: - if os.path.exists(drj_filename) and overwrite: - print('Removing and recreating {}'.format(drj_filename)) - os.remove(drj_filename) - science_hdu.writeto(drj_filename) - print('Defringed science data saved to {}'.format(drj_filename)) + if os.path.exists(output_filename) and overwrite: + print('Removing and recreating {}'.format(output_filename)) + os.remove(output_filename) + science_hdu.writeto(output_filename) + print('Defringed science data saved to {}'.format(output_filename)) - return drj_filename + return output_filename def parse_args(): @@ -190,7 +200,7 @@ def parse_args(): Output ---------------------------------------------------------------- Outputs a file that has the same rootname as the input science file, - but with '_drj.fits' at the end. This is the final defringed data. + but with '_defringe.fits' at the end. This is the final defringed data. '''), description='Script to calibrate science data in preparation for defringing') parser.add_argument('science', From 8da2ab6752ddb98bcc325c898b1eb6a2b3eebae3 Mon Sep 17 00:00:00 2001 From: Philip Hodge Date: Fri, 1 Nov 2019 15:27:13 -0400 Subject: [PATCH 6/6] Make the string for the regex a raw string. --- stistools/defringe.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stistools/defringe.py b/stistools/defringe.py index 034fe44a..5f856a67 100755 --- a/stistools/defringe.py +++ b/stistools/defringe.py @@ -62,7 +62,7 @@ def defringe(science_file, fringe_flat, output=None, science_file = os.path.normpath(expandFileName(science_file)) # Expand IRAF and UNIX $VARS sci_dir, sci_filename = os.path.split(science_file) if output is None: - sci_root = re.split('\.fits.*', sci_filename, + sci_root = re.split(r'\.fits.*', sci_filename, flags=re.IGNORECASE)[0].rsplit('_',1)[0] output_filename = os.path.join(sci_dir, sci_root + '_defringe.fits') else: