Check for NaN and comparison stars and remove them

stellarphot/differential_photometry/aij_rel_fluxes.py

@@ -27,7 +27,15 @@ def calc_aij_relative_flux(star_data, comp_stars,
         Table of star data from one or more images.
 
     comp_stars : '~astropy.table.Table'
-        Table of known comparison stars in the field, given by AAVSO
+        Table of comparison stars in the field. Must contain a column
+        called ``RA`` and a column called ``Dec``.
+        NOTE that not all
+        of the comparison stars will necessarily be used. Stars in
+        this table are excluded from the comparison set if, in any
+        of the `star_data` for that comparison, the net counts are
+        ``NaN`` or if the angular distance between the position in
+        the `star_data` and the position in the `comp_stars` table
+        is too large.
 
     in_place : ``bool``, optional
         If ``True``, add new columns to input table. Otherwise, return
@@ -90,8 +98,21 @@ def calc_aij_relative_flux(star_data, comp_stars,
     check_for_bad = check_for_bad.group_by('star_id')
     is_all_good = check_for_bad.groups.aggregate(np.all)
 
-    bad_comps = is_all_good['star_id'][~is_all_good['good']]
+    bad_comps = set(is_all_good['star_id'][~is_all_good['good']])
+    print(f'{bad_comps=}')
 
+    # Check whether any of the comp stars have NaN values and,
+    # if they do, exclude them from the comp set.
+    check_for_nan = Table(data=[star_data[star_id_column].data,
+                                star_data[flux_column_name].data],
+                          names=['star_id', 'net_counts'])
+    check_for_nan = check_for_nan.group_by('star_id')
+    check_for_nan['good'] = ~np.isnan(check_for_nan['net_counts'])
+    is_all_good = check_for_nan.groups.aggregate(np.all)
+
+    bad_comps = bad_comps | set(is_all_good['star_id'][~is_all_good['good']])
+
+    print(f'{bad_comps=}')
     for comp in bad_comps:
         this_comp = star_data[star_id_column] == comp
         good[this_comp] = False

stellarphot/differential_photometry/tests/test_aij_rel_fluxes.py

@@ -91,18 +91,27 @@ def test_relative_flux_calculation(in_place,
         assert 'relative_flux' not in input_table.colnames
 
 
-def test_mislocated_comp_star():
+@pytest.mark.parametrize('bad_thing', ['RA', 'NaN'])
+def test_bad_comp_star(bad_thing):
     expected_flux, expected_error, input_table, comp_star = \
         _raw_photometry_table()
-    # "Jiggle" one of the stars by moving it by a few arcsec in one image.
-    # We'll do it in the last image.
+    # We'll do modify the "bad" property for the last star in the last
+    # image.
 
     # First, let's sort so the row we want to modify is the last one
     input_table.sort(['date-obs', 'star_id'])
-    last_one = input_table[-1]
-    coord_inp = SkyCoord(ra=last_one['RA'], dec=last_one['Dec'], unit=u.degree)
-    coord_bad_ra = coord_inp.ra + 3 * u.arcsecond
-    input_table['RA'][-1] = coord_bad_ra.degree
+
+    # Force a copy of this row so we have access to the original values
+    last_one = Table(input_table[-1])
+
+    if bad_thing == 'RA':
+        # "Jiggle" one of the stars by moving it by a few arcsec in one image.
+        coord_inp = SkyCoord(ra=last_one['RA'], dec=last_one['Dec'],
+                             unit=u.degree)
+        coord_bad_ra = coord_inp.ra + 3 * u.arcsecond
+        input_table['RA'][-1] = coord_bad_ra.degree
+    elif bad_thing == 'NaN':
+        input_table['aperture_net_flux'][-1] = np.nan
 
     output_table = calc_aij_relative_flux(input_table, comp_star,
                                           in_place=False)
@@ -118,6 +127,9 @@ def test_mislocated_comp_star():
                             comp_star['aperture_net_flux'][1])
     new_expected_flux[2] = (comp_star['aperture_net_flux'][1] /
                             comp_star['aperture_net_flux'][0])
+
     new_expected_flux[3] = expected_flux[3]
+    if bad_thing == 'NaN':
+        new_expected_flux[3] = np.nan
 
     np.testing.assert_allclose(new_expected_flux, output_table['relative_flux'][-4:])