Coverage for src / lstautorta / Auto_Check_DL2.py: 0%

1#!/usr/bin/env python3 

2import argparse 

3import os 

4from os import listdir 

5from os.path import isfile, join 

6 

7import matplotlib.pyplot as plt 

8import numpy as np 

9from astropy.table import vstack 

10from ctapipe.io import read_table 

11from matplotlib.backends.backend_pdf import PdfPages 

12 

13parser = argparse.ArgumentParser( 

14 description="Automatic Script for the DL2 check", formatter_class=argparse.ArgumentDefaultsHelpFormatter 

15) 

16parser.add_argument("-d", "--directory", default="/fefs/onsite/pipeline/rta/data/", help="Directory for data") 

17parser.add_argument("-da", "--date", default="20230705", help="Date of the run to check") 

18parser.add_argument("-r", "--run-id", default="13600", help="run id to check") 

19# parser.add_argument("-add", "--add-string", default="reco/", help="add a string to the path") 

20parser.add_argument("-add", "--add-string", default="", help="add a string to the path") 

21 

22args = parser.parse_args() 

23config = vars(args) 

24 

25 

26def plot_variable(table, pdf, name): 

27 print(name) 

28 mask = table["is_good_event"] == 1 

29 table = table[mask] 

30 plt.hist(table[name], bins=100) 

31 plt.xlabel(name) 

32 pdf.savefig() 

33 plt.close() 

34 # plt.show() 

35 

36 

37def main(): 

38 mypath_run_dir = config["directory"] + config["date"] + "/" + config["run_id"] + "/" + config["add_string"] 

39 if not os.path.exists(mypath_run_dir + "plots"): 

40 os.mkdir(mypath_run_dir + "plots") 

41 mypath = mypath_run_dir + "DL2/" 

42 filename = [f for f in listdir(mypath) if isfile(join(mypath, f)) if "dl2_v06_" in f] 

43 params = [] 

44 for i in range(len(filename)): 

45 tablename = "/dl2/event/telescope/parameters/LST_LSTCam" 

46 params.append(read_table(mypath + filename[i], tablename)) 

47 params = vstack(params) 

48 max_trigger_time = 1894367064 

49 mask = params["trigger_time"] < max_trigger_time # *(params['trigger_time']>1651613520) 

50 

51 start_time = np.sort(params["trigger_time"][mask])[0] 

52 end_time = np.sort(params["trigger_time"][mask])[-1] 

53 times = end_time - start_time 

54 print(times) 

55 nbins = int(times) 

56 filter_ = np.diff(params["event_id"][mask]) < 10000 

57 with PdfPages(mypath_run_dir + "plots/output_DL2.pdf") as pdf: 

58 plot_variable(params, pdf, "log_reco_energy") 

59 plot_variable(params, pdf, "reco_energy") 

60 plot_variable(params, pdf, "reco_disp_dx") 

61 plot_variable(params, pdf, "reco_disp_dy") 

62 plot_variable(params, pdf, "reco_src_x") 

63 plot_variable(params, pdf, "reco_src_y") 

64 plot_variable(params, pdf, "gammaness") 

65 

66 mask = params["is_good_event"] == 1 

67 params = params[mask] 

68 mask = params["gammaness"] > 0.7 

69 

70 plt.hist2d(params["reco_src_x"][mask], params["reco_src_y"][mask], bins=100) 

71 plt.xlabel("reco_src_x") 

72 plt.ylabel("reco_src_y") 

73 pdf.savefig() 

74 plt.close() 

75 

76 

77if __name__ == "__main__": 

78 main()