Coverage for lst_auto_rta/Plot_event_DL1.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 

9import tables 

10from matplotlib.backends.backend_pdf import PdfPages 

11 

12parser = argparse.ArgumentParser( 

13 description="Automatic Script for the DL1 check", formatter_class=argparse.ArgumentDefaultsHelpFormatter 

14) 

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

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

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

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

19parser.add_argument("-o", "--output", default="", help="output path") 

20 

21args = parser.parse_args() 

22config = vars(args) 

23 

24font = {"family": "normal", "weight": "bold", "size": 22} 

25 

26plt.rc("font", **font) 

27 

28 

29def plotImage(outputPlot, tabX, tabY, signal, imageIndex, eventId, eventType, tabFixPixelOrder=None, id_=0): 

30 strInfo = f"image {imageIndex}, id = {eventId}, type = {eventType},\n Signal mean = {signal.mean()} pe, min = {signal.min()}, max = {signal.max()}, std = {signal.std()}" 

31 print(strInfo) 

32 fig = plt.figure(figsize=(16, 10)) 

33 fig.patch.set_alpha(1.0) 

34 #print(len(tabFixPixelOrder)) 

35 

36 if tabFixPixelOrder is not None: 

37 plt.scatter(tabX[tabFixPixelOrder], tabY[tabFixPixelOrder], c=signal[:-7], s=120) 

38 else: 

39 plt.scatter(tabX, tabY, c=signal, s=120) 

40 

41 plt.axis("equal") 

42 plt.xlabel("x") 

43 plt.ylabel("y") 

44 plt.colorbar() 

45 plt.text(0.05, 0.95, strInfo, transform=fig.transFigure, size=14) 

46 plt.show() 

47 plt.savefig(outputPlot + "/plotImage" + str(id_) + ".png") 

48 # pdf.savefig() # saves the current figure into a pdf page 

49 plt.close() 

50 

51 

52def processAllFile(outputPlot, inputFile, firstImageIndex=0, lastImageIndex=10000, indexStep=1, id_=0): 

53 hfile = tables.open_file(inputFile, "r") 

54 

55 tabX = hfile.root.configuration.instrument.telescope.camera.geometry_LSTCam.col("pix_x") 

56 tabY = hfile.root.configuration.instrument.telescope.camera.geometry_LSTCam.col("pix_y") 

57 

58 try: 

59 tabPixelOrder = hfile.root.configuration.instrument.telescope.camera.pixel_order.tel_001.read() 

60 except: 

61 tabPixelOrder = None 

62 #print(hfile.root.r0.monitoring.telescope.unusable_pixels.tel_001[0]) 

63 # tabImages = hfile.root.dl1.event.telescope.images.tel_001.col("image") 

64 # tabEventId = hfile.root.dl1.event.subarray.trigger.col("event_id") 

65 # tabEventType = hfile.root.dl1.event.subarray.trigger.col("event_type") 

66 

67 tabImages = hfile.root.dl1.event.telescope.images.tel_001.col("image") 

68 tabTime = hfile.root.dl1.event.telescope.images.tel_001.col("peak_time") 

69 tabEventId = hfile.root.dl1.event.telescope.parameters.tel_001.col("event_id") 

70 tabEventType = hfile.root.dl1.event.telescope.parameters.tel_001.col("event_type") 

71 tabIntensity = hfile.root.dl1.event.telescope.parameters.tel_001.col("hillas_intensity") 

72 # id_ = 0 

73 # for i in range(firstImageIndex, len(tabImages), indexStep): 

74 for i in range(firstImageIndex, lastImageIndex, indexStep): 

75 if tabEventType[i] == 32 and tabImages[i].std() > 5: #and tabIntensity[i] < 12000: 

76 plotImage(outputPlot, tabX, tabY, tabImages[i], i, tabEventId[i], tabEventType[i], tabPixelOrder, id_) 

77 id_ = id_ + 1 

78 return id_ 

79 

80 

81mypath = config["directory"] + config["date"] + "/" + config["run_id"] + "/" + config["add_string"] + "/DL1/" 

82mypath_rundir = config["directory"] + config["date"] + "/" + config["run_id"] + "/" + config["add_string"] 

83print("Path is here ", mypath) 

84if not os.path.exists(mypath_rundir + "plots"): 

85 os.mkdir(mypath_rundir + "plots") 

86filename = [f for f in listdir(mypath) if isfile(join(mypath, f)) if "dl1" in f if not "v06" in f] 

87print("Plotting image from these files", filename) 

88index = 0 

89for i in filename: 

90 print(mypath + i) 

91 print("output in ", mypath_rundir + "plots") 

92 index = processAllFile(mypath_rundir + "plots", mypath + i, id_=index)