Skip to content

Commit

Permalink
Merge pull request #468 from mwcraig/revive-transforms
Browse files Browse the repository at this point in the history
Revive transforms
  • Loading branch information
mwcraig authored Oct 10, 2024
2 parents cb1dcda + 5260c1c commit adc4858
Show file tree
Hide file tree
Showing 3 changed files with 155 additions and 28 deletions.
121 changes: 100 additions & 21 deletions stellarphot/notebooks/photometry/04-transform-pared-back.ipynb
100755 → 100644
Original file line number Diff line number Diff line change
Expand Up @@ -19,11 +19,13 @@
"import matplotlib.pyplot as plt\n",
"\n",
"from stellarphot.utils.magnitude_transforms import (\n",
" f,\n",
" calibrated_from_instrumental,\n",
" opts_to_str,\n",
" calc_residual,\n",
" transform_to_catalog,\n",
")"
")\n",
"\n",
"from stellarphot import PhotometryData, vsx_vizier"
]
},
{
Expand Down Expand Up @@ -74,29 +76,49 @@
"c_min = -0.5\n",
"d_min = -1e-6\n",
"\n",
"our_filters = [\"SI\"]\n",
"our_filters = [\"SR\"]\n",
"\n",
"aavso_band_names = dict(B=\"B\", V=\"V\", gp=\"SG\", rp=\"SR\", ip=\"SI\", SI=\"SI\", SR=\"SR\")\n",
"\n",
"aavso_band_names = dict(B=\"B\", V=\"V\", gp=\"SG\", rp=\"SR\", ip=\"SI\", SI=\"SI\")\n",
"# cat_color_colums = dict(\n",
"# B=(\"Bmag\", \"Vmag\"),\n",
"# V=(\"Bmag\", \"Vmag\"),\n",
"# gp=(\"g_mag\", \"r_mag\"),\n",
"# rp=(\"r_mag\", \"i_mag\"),\n",
"# ip=(\"r_mag\", \"i_mag\"),\n",
"# SI=(\"r_mag\", \"i_mag\"),\n",
"# )\n",
"\n",
"cat_color_colums = dict(\n",
" B=(\"Bmag\", \"Vmag\"),\n",
" V=(\"Bmag\", \"Vmag\"),\n",
" gp=(\"g_mag\", \"r_mag\"),\n",
" rp=(\"r_mag\", \"i_mag\"),\n",
" ip=(\"r_mag\", \"i_mag\"),\n",
" SI=(\"r_mag\", \"i_mag\"),\n",
" B=(\"mag_B\", \"mag_V\"),\n",
" V=(\"mag_B\", \"mag_V\"),\n",
" gp=(\"mag_SG\", \"mag_SR\"),\n",
" rp=(\"mag_SR\", \"mag_SI\"),\n",
" SR=(\"SR\", \"SI\"),\n",
" ip=(\"mag_SR\", \"mag_SI\"),\n",
" SI=(\"mag_SR\", \"mag_SI\"),\n",
")\n",
"\n",
"# cat_filter = dict(\n",
"# B=\"Bmag\",\n",
"# V=\"Vmag\",\n",
"# gp=\"g_mag\",\n",
"# rp=\"r_mag\",\n",
"# ip=\"i_mag\",\n",
"# SI=\"i_mag\",\n",
"# )\n",
"\n",
"cat_filter = dict(\n",
" B=\"Bmag\",\n",
" V=\"Vmag\",\n",
" gp=\"g_mag\",\n",
" rp=\"r_mag\",\n",
" ip=\"i_mag\",\n",
" SI=\"i_mag\",\n",
" B=\"B\",\n",
" V=\"mag_V\",\n",
" gp=\"mag_SG\",\n",
" rp=\"mag_SR\",\n",
" SR=\"SR\",\n",
" ip=\"mag_SI\",\n",
" SI=\"mag_SI\",\n",
")\n",
"\n",
"input_photometry_file = \"TIC-402828941-2021-09-28.ecsv\"\n",
"input_photometry_file = \"sp2-kelt-1b-photometry-2021-09-14.ecsv\"\n",
"output_photometry_file = \"some_name.ecsv\""
]
},
Expand All @@ -108,7 +130,11 @@
},
"outputs": [],
"source": [
"all_mags = Table.read(input_photometry_file)\n",
"all_mags = PhotometryData.read(input_photometry_file)\n",
"\n",
"# BAD BAD BAD BAD BAD BAD BAD BAD BAD BAD\n",
"all_mags[\"passband\"] = \"SR\"\n",
"\n",
"\n",
"# Ensure we have the right table ordering later\n",
"all_mags.sort([\"passband\", \"bjd\"])"
Expand Down Expand Up @@ -137,6 +163,15 @@
"assert all(k[0] in cat_filter.keys() for k in filter_groups.groups.keys)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"filter_groups.groups.keys"
]
},
{
"cell_type": "markdown",
"metadata": {},
Expand All @@ -148,6 +183,7 @@
"cell_type": "code",
"execution_count": null,
"metadata": {
"scrolled": true,
"tags": []
},
"outputs": [],
Expand All @@ -156,7 +192,7 @@
"\n",
"for k, group in zip(our_filters, filter_groups.groups):\n",
" print(f\"Transforming band {k}\")\n",
" by_bjd = group.group_by(\"bjd\")\n",
" by_bjd = group.group_by(\"file\")\n",
"\n",
" transform_to_catalog(\n",
" by_bjd,\n",
Expand Down Expand Up @@ -190,7 +226,32 @@
"metadata": {},
"outputs": [],
"source": [
"plt.plot(by_bjd[\"mag_cat\"] - by_bjd[\"mag_inst_cal\"], \"o\")\n",
"plt.plot(by_bjd[\"mag_cat\"], by_bjd[\"mag_cat\"] - by_bjd[\"mag_inst_cal\"], \"o\")\n",
"plt.ylim(-0.05, 0.05)\n",
"plt.grid()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plt.figure()\n",
"kelt1 = by_bjd.lightcurve_for(\"kelt-1\", flux_column=\"mag_inst_cal\") # by_bjd[\"star_id\"] == 1\n",
"plt.plot(kelt1.time.jd, kelt1.flux, \".\")\n",
"plt.ylim(reversed(plt.ylim()))\n",
"plt.grid()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plt.figure()\n",
"plt.plot(kelt1.time.value, kelt1[\"z\"], \".\")\n",
"plt.grid()"
]
},
Expand All @@ -203,6 +264,24 @@
"output_table.write(output_photometry_file)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"by_bjd.colnames"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"by_bjd['id']"
]
},
{
"cell_type": "code",
"execution_count": null,
Expand All @@ -227,7 +306,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.5"
"version": "3.11.9"
}
},
"nbformat": 4,
Expand Down
58 changes: 53 additions & 5 deletions stellarphot/utils/magnitude_transforms.py
Original file line number Diff line number Diff line change
Expand Up @@ -424,7 +424,7 @@ def transform_magnitudes(
transform_catalog_coords = SkyCoord(
transform_catalog["RAJ2000"], transform_catalog["DEJ2000"], unit="deg"
)
input_coords = SkyCoord(input_mags["RA"], input_mags["Dec"])
input_coords = SkyCoord(input_mags["ra"], input_mags["dec"])

transform_catalog_index, d2d, _ = match_coordinates_sky(
input_coords, transform_catalog_coords
Expand Down Expand Up @@ -593,22 +593,51 @@ def transform_to_catalog(
cat = None
cat_coords = None

def fake_it(one_image):
# Accumulate the parameters
popt = [np.nan] * len(all_params)
for param, value in zip(all_params, popt, strict=True):
param.extend([value] * len(one_image))
nana = [np.nan] * len(one_image)
cal_mags.extend(nana)
cat_mags.extend(nana)
cat_colors.extend(nana)
resids.extend(nana)

for file, one_image_inp in zip(
observed_mags_grouped.groups.keys, observed_mags_grouped.groups, strict=True
):
one_image = one_image_inp[one_image_inp["passband"] == obs_filter]
our_coords = SkyCoord(one_image["RA"], one_image["Dec"], unit="degree")
our_coords = SkyCoord(one_image["ra"], one_image["dec"], unit="degree")
if cat is None or cat_coords is None:
cat = apass_dr9(
our_coords[0], radius=1 * u.degree, clip_by_frame=False, padding=0
)
cat_coords = SkyCoord(cat["ra"], cat["dec"])
cat["color"] = cat[cat_color[0]] - cat[cat_color[1]]
cat = cat.passband_columns(passbands=["B", "V", "SR", "SG", "SI"])
cat["mag_RC"] = filter_transform(
cat,
output_filter="R",
g="mag_SG",
r="mag_SR",
i="mag_SI",
transform="jester",
)
cat["mag_IC"] = filter_transform(
cat,
output_filter="I",
g="mag_SG",
r="mag_SR",
i="mag_SI",
transform="jester",
)

cat_coords = SkyCoord(cat["ra"], cat["dec"], unit="degree")
cat["color"] = cat[f"mag_{cat_color[0]}"] - cat[f"mag_{cat_color[1]}"]

cat_idx, d2d, _ = our_coords.match_to_catalog_sky(cat_coords)

mag_inst = one_image[obs_mag_col]
cat_mag = cat[cat_filter][cat_idx]
cat_mag = cat[f"mag_{cat_filter}"][cat_idx]
color = cat["color"][cat_idx]

# Impose some constraints on what is included in the fit
Expand All @@ -619,8 +648,22 @@ def transform_to_catalog(
& ~np.isnan(one_image[obs_mag_col])
)

if not (any(good_dat) or any(good_cat)):
print(f"No good data in {file[0]}")
fake_it(one_image)
continue

mag_diff = cat_mag - mag_inst

if (
(not (any(good_dat) or any(good_cat)))
or all(np.isnan(mag_diff))
or all(mag_diff.mask)
):
print(f"No good data in {file[0]}")
fake_it(one_image)
continue

good_data = good_dat & (
np.abs(mag_diff - np.nanmedian(mag_diff[good_dat & ~mag_diff.mask])) < 1
)
Expand All @@ -631,6 +674,11 @@ def transform_to_catalog(

good = good_cat & good_data

if not any(good):
print(f"No good data in {file[0]}")
fake_it(one_image)
continue

# Prep for fitting
init_guess = (a_cen, 0, 0, 0, zero_point_mid)
X = (mag_inst[good], color[good])
Expand Down
4 changes: 2 additions & 2 deletions stellarphot/utils/tests/test_magnitude_transforms.py
Original file line number Diff line number Diff line change
Expand Up @@ -73,8 +73,8 @@ def generate_tables(n_stars, mag_model):
ra, dec = generate_star_coordinates(n_stars)

# Instrumental magnitudes
ra_col = Column(name="RA", data=ra)
dec_col = Column(name="Dec", data=dec)
ra_col = Column(name="ra", data=ra)
dec_col = Column(name="dec", data=dec)

instrumental = Table([instr_mags, ra_col, dec_col])

Expand Down

0 comments on commit adc4858

Please sign in to comment.