From 39417633b19016ffc8ff21154f2e7b51eb97d115 Mon Sep 17 00:00:00 2001 From: Elisa-Visentin Date: Thu, 12 Sep 2024 13:11:43 +0000 Subject: [PATCH 1/3] ra dec db --- .../database_production/create_LST_table.py | 8 + .../database_production/set_ra_dec.py | 171 ++++++++++++++++++ setup.cfg | 1 + 3 files changed, 180 insertions(+) create mode 100644 magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py diff --git a/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/create_LST_table.py b/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/create_LST_table.py index 9fd1c0696..0d24c361e 100644 --- a/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/create_LST_table.py +++ b/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/create_LST_table.py @@ -88,6 +88,14 @@ def main(): out_h5, key=out_key, ) + if "ra" in df_old: + df_cut["ra"] = np.nan + if "dec" in df_old: + df_cut["dec"] = np.nan + if "MC_dec" in df_old: + df_cut["MC_dec"] = np.nan + if "point_source" in df_old: + df_cut["point_source"] = np.nan df_cut = pd.concat([df_old, df_cut]).drop_duplicates( subset="LST1_run", keep="first" ) diff --git a/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py b/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py new file mode 100644 index 000000000..61425e5f4 --- /dev/null +++ b/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py @@ -0,0 +1,171 @@ +""" +Add, to LST database, infos about coordinates and extension of the sources and MC declination to be used to process the source + +Usage: +$ set_ra_dec (-b YYYYMMDD -e YYYYMMDD -r ra_dict -d dec_dict -m mc_dec -p point_source_dict) +""" + +import argparse +import json + +import numpy as np +import pandas as pd +import yaml +from astropy.coordinates import SkyCoord +from astropy.coordinates.name_resolve import NameResolveError + +from magicctapipe.io import resource_file + + +def main(): + + """ + Main function + """ + + parser = argparse.ArgumentParser() + + parser.add_argument( + "--begin-date", + "-b", + dest="begin", + type=int, + default=0, + help="First date to update database (YYYYMMDD)", + ) + parser.add_argument( + "--end-date", + "-e", + dest="end", + type=int, + default=0, + help="End date to update database (YYYYMMDD)", + ) + parser.add_argument( + "--dict-source", + "-s", + dest="source", + type=str, + default="./source_dict.json", + help="File with dictionary of info (RA/Dec/point-source) for sources", + ) + parser.add_argument( + "--mc-dec", + "-m", + dest="dec_mc", + type=str, + default="./dec_mc.json", + help="File with list of MC declinations", + ) + + args = parser.parse_args() + config_file = resource_file("database_config.yaml") + + with open( + config_file, "rb" + ) as fc: # "rb" mode opens the file in binary format for reading + config_dict = yaml.safe_load(fc) + + LST_h5 = config_dict["database_paths"]["LST"] + LST_key = config_dict["database_keys"]["LST"] + + df_LST = pd.read_hdf( + LST_h5, + key=LST_key, + ) + if "ra" not in df_LST: + df_LST["ra"] = np.nan + if "dec" not in df_LST: + df_LST["dec"] = np.nan + if "MC_dec" not in df_LST: + df_LST["MC_dec"] = np.nan + if "point_source" not in df_LST: + df_LST["point_source"] = np.nan + df_LST_full = df_LST.copy(deep=True) + if args.begin != 0: + df_LST = df_LST[df_LST["DATE"].astype(int) >= args.begin] + if args.end != 0: + df_LST = df_LST[df_LST["DATE"].astype(int) <= args.end] + + sources = np.unique(df_LST["source"]) + with open(args.source) as f: + dict_source = f.read() + with open(args.dec_mc) as f: + mc_dec = f.read() + source_dict = json.loads(dict_source) + dec_mc = np.asarray(json.loads(mc_dec)).astype(np.float64) + print("MC declinations: \t", dec_mc) + print("\n\nChecking RA/Dec...\n\n") + i = 0 + for src in sources: + + try: + coord = SkyCoord.from_name(src) + if src == "Crab": + coord = SkyCoord.from_name("CrabNebula") + src_dec = coord.dec.degree + src_ra = coord.ra.degree + df_LST["ra"] = np.where(df_LST["source"] == src, src_ra, df_LST["ra"]) + df_LST["dec"] = np.where(df_LST["source"] == src, src_dec, df_LST["dec"]) + df_LST["MC_dec"] = np.where( + df_LST["source"] == src, + float(dec_mc[np.argmin(np.abs(src_dec - dec_mc))]), + df_LST["MC_dec"], + ) + except NameResolveError: + print(f"{i}: {src} not found in astropy. Looking to the dictionaries...") + if (src in source_dict.keys()) and (source_dict.get(src)[0] != "NaN"): + src_ra = float(source_dict.get(src)[0]) + df_LST["ra"] = np.where(df_LST["source"] == src, src_ra, df_LST["ra"]) + + else: + print( + f"\t {i}: {src} RA not in the dictionaries. Please add it to the dictionaries" + ) + + if (src in source_dict.keys()) and (source_dict.get(src)[1] != "NaN"): + src_dec = float(source_dict.get(src)[1]) + df_LST["dec"] = np.where( + df_LST["source"] == src, src_dec, df_LST["dec"] + ) + df_LST["MC_dec"] = np.where( + df_LST["source"] == src, + float(dec_mc[np.argmin(np.abs(src_dec - dec_mc))]), + df_LST["MC_dec"], + ) + else: + print( + f"\t {i}: {src} Dec not in the dictionaries. Please add it to the dictionaries" + ) + + i += 1 + print("\n\nChecking if point source...\n\n") + i = 0 + for src in sources: + if (src in source_dict.keys()) and (source_dict.get(src)[2] != "NaN"): + src_point = str(source_dict.get(src)[2]) + df_LST["point_source"] = np.where( + df_LST["source"] == src, src_point, df_LST["point_source"] + ) + else: + print( + f"\t {i}: {src} extension information not in the dictionaries. Please add it to the dictionaries" + ) + i += 1 + df_LST = pd.concat([df_LST, df_LST_full]).drop_duplicates( + subset="LST1_run", keep="first" + ) + df_LST.to_hdf( + LST_h5, + key=LST_key, + mode="w", + min_itemsize={ + "lstchain_versions": 20, + "last_lstchain_file": 90, + "processed_lstchain_file": 90, + }, + ) + + +if __name__ == "__main__": + main() \ No newline at end of file diff --git a/setup.cfg b/setup.cfg index ecd18a554..5ec0dfeea 100644 --- a/setup.cfg +++ b/setup.cfg @@ -108,6 +108,7 @@ console_scripts = merging_runs = magicctapipe.scripts.lst1_magic.semi_automatic_scripts.merging_runs:main nsb_level = magicctapipe.scripts.lst1_magic.semi_automatic_scripts.database_production.nsb_level:main nsb_to_h5 = magicctapipe.scripts.lst1_magic.semi_automatic_scripts.database_production.nsb_to_h5:main + set_ra_dec = magicctapipe.scripts.lst1_magic.semi_automatic_scripts.database_production.set_ra_dec:main stereo_events = magicctapipe.scripts.lst1_magic.semi_automatic_scripts.stereo_events:main update_MAGIC_database = magicctapipe.scripts.lst1_magic.semi_automatic_scripts.database_production.update_MAGIC_database:main From 10c55ea33c7c101e24551c30a9d80bb806dd7ef9 Mon Sep 17 00:00:00 2001 From: Elisa-Visentin Date: Thu, 12 Sep 2024 13:26:36 +0000 Subject: [PATCH 2/3] minor fixes --- .../semi_automatic_scripts/database_production/set_ra_dec.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py b/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py index 61425e5f4..654692148 100644 --- a/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py +++ b/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py @@ -2,7 +2,7 @@ Add, to LST database, infos about coordinates and extension of the sources and MC declination to be used to process the source Usage: -$ set_ra_dec (-b YYYYMMDD -e YYYYMMDD -r ra_dict -d dec_dict -m mc_dec -p point_source_dict) +$ set_ra_dec (-b YYYYMMDD -e YYYYMMDD -s source_dict -m mc_dec) """ import argparse @@ -168,4 +168,4 @@ def main(): if __name__ == "__main__": - main() \ No newline at end of file + main() From 027529b56e7f5c9871208f85bd3f56dd25c82e42 Mon Sep 17 00:00:00 2001 From: Elisa-Visentin Date: Fri, 13 Sep 2024 09:49:37 +0000 Subject: [PATCH 3/3] minor fixes --- .../database_production/set_ra_dec.py | 43 ++++++++----------- 1 file changed, 18 insertions(+), 25 deletions(-) diff --git a/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py b/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py index 654692148..ac2595687 100644 --- a/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py +++ b/magicctapipe/scripts/lst1_magic/semi_automatic_scripts/database_production/set_ra_dec.py @@ -105,40 +105,33 @@ def main(): coord = SkyCoord.from_name("CrabNebula") src_dec = coord.dec.degree src_ra = coord.ra.degree - df_LST["ra"] = np.where(df_LST["source"] == src, src_ra, df_LST["ra"]) - df_LST["dec"] = np.where(df_LST["source"] == src, src_dec, df_LST["dec"]) - df_LST["MC_dec"] = np.where( - df_LST["source"] == src, - float(dec_mc[np.argmin(np.abs(src_dec - dec_mc))]), - df_LST["MC_dec"], - ) + except NameResolveError: print(f"{i}: {src} not found in astropy. Looking to the dictionaries...") - if (src in source_dict.keys()) and (source_dict.get(src)[0] != "NaN"): + if ( + (src in source_dict.keys()) + and (source_dict.get(src)[0] != "NaN") + and (source_dict.get(src)[1] != "NaN") + ): src_ra = float(source_dict.get(src)[0]) - df_LST["ra"] = np.where(df_LST["source"] == src, src_ra, df_LST["ra"]) - - else: - print( - f"\t {i}: {src} RA not in the dictionaries. Please add it to the dictionaries" - ) - - if (src in source_dict.keys()) and (source_dict.get(src)[1] != "NaN"): src_dec = float(source_dict.get(src)[1]) - df_LST["dec"] = np.where( - df_LST["source"] == src, src_dec, df_LST["dec"] - ) - df_LST["MC_dec"] = np.where( - df_LST["source"] == src, - float(dec_mc[np.argmin(np.abs(src_dec - dec_mc))]), - df_LST["MC_dec"], - ) + else: print( - f"\t {i}: {src} Dec not in the dictionaries. Please add it to the dictionaries" + f"\t {i}: {src} RA and/or Dec not in the dictionary. Please update the dictionary" ) + src_ra = np.nan + src_dec = np.nan i += 1 + df_LST["ra"] = np.where(df_LST["source"] == src, src_ra, df_LST["ra"]) + df_LST["dec"] = np.where(df_LST["source"] == src, src_dec, df_LST["dec"]) + if not (np.isnan(src_dec)): + df_LST["MC_dec"] = np.where( + df_LST["source"] == src, + float(dec_mc[np.argmin(np.abs(src_dec - dec_mc))]), + df_LST["MC_dec"], + ) print("\n\nChecking if point source...\n\n") i = 0 for src in sources: