import argparse import os import matplotlib.pyplot as plt import numpy as np from glob import glob import matplotlib as mpl import matplotlib.gridspec as gridspec import re from readers.topreader import TopReader from plotters.histoplotter import HistoPlotter from typing import Dict from readers.mocanloreader import MocanloReader from copy import deepcopy ########################################## ####### plot names and setup ############# ########################################## with open("MoCaNLO_STXS_binning.txt", 'r') as file: lines = file.readlines() # Reads all lines into a list merge_DPhijj_pTHjj=False # read the text file of the STXS binning STXS_bin_names = [] for line in lines: if line=="\n" or "#" in line: continue new_line = line.strip("\n") new_line = new_line.split(".")[1].strip() STXS_bin_names.append(new_line) keys_to_plot = {"histogram_Mjj_pTH_truth_Mjj_vs_pTH": ["MJJ-PTH-0-80-PTJ-","MJJ-PTH-80-120-PTJ-","MJJ-PTH-120-260-PTJ-", "MJJ-PTH-260-500-PTJ-","MJJ-PTH-500-850-PTJ-","MJJ-PTH-850-INFTY-PTJ-"], "histogram_Mjj_Dphijj_truth_Mjj_vs_Dphijj":["MJJ-DPHIJJ-0-PIov4-PTJ-", "MJJ-DPHIJJ-PIov4-PIov2-PTJ-", "MJJ-DPHIJJ-PIov2-3PIov4-PTJ-", "MJJ-DPHIJJ-3PIov4-PI-PTJ-"], "histogram_Mjj_Dyjj_truth_Mjj_vs_Dyjj": ["MJJ-DYJJ-2-4-PTJ-", "MJJ-DYJJ-4-5-PTJ-","MJJ-DYJJ-5-6-PTJ-","MJJ-DYJJ-6-7-PTJ-","MJJ-DYJJ-7-INFTY-PTJ-"], "histogram_pTH_Dyjj_truth_pTH_vs_Dyjj":["PTH-DYJJ-2-4-PTJ-","PTH-DYJJ-4-5-PTJ-","PTH-DYJJ-5-6-PTJ-","PTH-DYJJ-6-7-PTJ-","PTH-DYJJ-7-INFTY-PTJ-"], "histogram_Mjj_pTH_STSX_truth_Mjj_vs_pTH_STSX": STXS_bin_names, "histogram_Mjj_pTH_STSX_truth_Mjj_vs_pTH_STSX_Mjj_below_350": STXS_bin_names} fiducial_setups = {"20": "setup.2.a", "30": "setup.2.b" } stxs_setup = {"stxs": "STXS"} merge_pTHjj = True ########################################## ########### plotting settings ############# ########################################## proVBF_setup = "histo" default_label_mapping = { "full_nlo_ew_photon_only": "Photon only", "full_nlo_ew_w_photon": "NLO with photon", "best_prediction_final": "Best prediction", "MoCaNLO_NLO_EW_photon_only": "Photon only", "best_prediction_1": "without photon", "best_prediction_2" : "photon only correction", "best_prediction_3" : "EW NLO with photon - EW NLO without photon", "bigstxs_postprocessed_proVBFH_LO": "LO", "bigstxs_postprocessed_proVBFH_NLO": "NLO QCD", "STXS_proVBFH_NNLO_EW": "Best prediction", "bigstxs_postprocessed_proVBFH_LO_ratio": "LO", "bigstxs_postprocessed_proVBFH_NLO_ratio": "NLO", "bigstxs_postprocessed_proVBFH_NNLO_ratio": "NNLO", "bigstxs_postprocessed_proVBFH_NNLO": "NNLO QCD", "%s_central_scale_proVBFH_LO"%proVBF_setup: "LO", "%s_central_scale_proVBFH_NLO"%proVBF_setup: "NLO", "%s_central_scale_proVBFH_NNLO"%proVBF_setup: "NNLO", "%s_proVBFH_NNLO_nonfact_ratio"%proVBF_setup: "NF", "%s_central_scale_proVBFH_LO_ratio"%proVBF_setup: "LO", "%s_central_scale_proVBFH_NLO_ratio"%proVBF_setup: "NLO", "%s_central_scale_proVBFH_NNLO_ratio"%proVBF_setup: "NNLO", "full_nlo_qcd_correction": "NLO", "mocanlo_qcd_nlo_prediction": "NLO", "mocanlo_nlo_qcd_born": "LO", "MoCaNLO_NLO_EW": "", "MoCaNLO_NLO_QCD_ratio": "NLO", "MoCaNLO_LO_QCD_ratio": "LO", "MoCaNLO_loop_induced": "Loop-induced interference", "MoCaNLO_loop_induced_2": "Loop-induced squared", "MoCaNLO_ggF_LO": "LO ggF", "MoCaNLO_NLO_QCD_ratio_Mjj_smaller_350": "NLO", "MoCaNLO_LO_QCD_ratio_Mjj_smaller_350": "LO" } default_color_mapping = { "full_nlo_ew_photon_only": "#E1C16E", #Brass "full_nlo_ew_w_photon": "#6F4E37", #coffee "best_prediction_final": "#50C878", #Emerald Green "MoCaNLO_NLO_EW_photon_only": "", "best_prediction_1": "", "best_prediction_2" : "", "best_prediction_3" : "", "%s_central_scale_proVBFH_LO"%proVBF_setup: "#0096FF", #bright blue "%s_central_scale_proVBFH_NLO"%proVBF_setup: "#5F9EA0", #Cadet Blue "STXS_proVBFH_NNLO_EW": "", "%s_central_scale_proVBFH_LO_ratio"%proVBF_setup: "#0096FF", "%s_central_scale_proVBFH_NLO_ratio"%proVBF_setup: "#5F9EA0", #Cadet Blue "MoCaNLO_LO_QCD_ratio_Mjj_smaller_350": "#0096FF", "MoCaNLO_NLO_QCD_ratio_Mjj_smaller_350": "#5F9EA0", #Cadet Blue "full_nlo_qcd_correction": "#5F9EA0", #Cadet Blue "mocanlo_qcd_nlo_prediction": "#5F9EA0", "mocanlo_nlo_qcd_born": "#0096FF", "%s_central_scale_proVBFH_NNLO_ratio"%proVBF_setup: "#C30000", #"#5D3FD3",#Iris "%s_proVBFH_NNLO_nonfact_ratio"%proVBF_setup: "#F5CF9F", #light orange "MoCaNLO_NLO_EW": "", "MoCaNLO_NLO_QCD_ratio": "#988558", #Dark Tan "MoCaNLO_LO_QCD_ratio": "#CD7F32", #Bronze "MoCaNLO_loop_induced": "#708090",#Slate Gray "MoCaNLO_loop_induced_2": "#8A9A5B",#Sage Green "MoCaNLO_ggF_LO":"#F5CF9F" } ####specific setting for the different 2D or STXS plots # reference_list is just a settings json for plotting #default is for pTH vs Mjj default_reference_list = [ {"data_ind": ["%s_central_scale_proVBFH_LO"%proVBF_setup, "best_prediction_final"],#"best_prediction_1", "best_prediction_2", "best_prediction_3"], "error_bands": True, #"y_scale": (10**(-8),5*10**(-3)), }, {"data_ind": ["%s_central_scale_proVBFH_NLO_ratio"%proVBF_setup, "%s_central_scale_proVBFH_LO_ratio"%proVBF_setup, "%s_central_scale_proVBFH_NNLO_ratio"%proVBF_setup, "%s_proVBFH_NNLO_nonfact_ratio"%proVBF_setup], #"reference_ind": "bigstxs_postprocessed_proVBFH_NLO", "error_bands": True, "horizontal_line": False, #"plot_label_for_reference": True, "label": "ratio to\nVBF NLO QCD", "y_scale": (0.4,1.53), }, {"data_ind": ["full_nlo_ew_photon_only", "full_nlo_ew_w_photon"], "error_bands": False, "horizontal_line": True, "label": "Full\n NLO EW/LO", "y_scale": (0.65,1.02), }, {"data_ind": ["MoCaNLO_NLO_QCD_ratio", "MoCaNLO_LO_QCD_ratio"], "error_bands": True, "label": "Full/VBF", "horizontal_line": True, "y_scale": (0.0000,6.3000), }, {"data_ind": ["MoCaNLO_loop_induced", "MoCaNLO_loop_induced_2"],#, "MoCaNLO_ggF_LO"], "reference_ind": "MoCaNLO_LO", "horizontal_line": False, "label": "Loop-induced/LO", "y_scale": (-0.005,0.7000), "error_bands": True, }, {"data_ind": ["MoCaNLO_ggF_LO"],#, "MoCaNLO_ggF_LO"], "reference_ind": "MoCaNLO_LO", "horizontal_line": False, "label": "ggF/LO", "y_scale": (-0.0,10.5000), "error_bands": True, } ] small_ticks = ["[300,500)","[500,700)","[700,900)","[900,1100)","[1100,$\\infty$)"] medium_ticks = ["$p_T^H \\in$[0,80)","$p_T^H \\in$[80,120)","$p_T^H \\in$[120,260)","$p_T^H \\in$[260,500)","$p_T^H \\in$[500,850)","$p_T^H \\in$[850,$\\infty$)"] Nbins = 4 default_settings = { "complex_indexes": { "small_ticks": small_ticks, "small_ticks_rotation": 65, "medium_ticks": medium_ticks, "medium_ticks_labelpad": 90, "big_ticks_labelpad": 120, "big_ticks": [""], #"v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "big_ticks_xlines": [], }, "x_param_labelpad": 60 } default_figsize = (14,11) def plotting_setup(key_to_plot, proVBF_setup, Nbins): reference_list = [ {"data_ind": ["%s_central_scale_proVBFH_LO"%proVBF_setup, "best_prediction_final"],#"best_prediction_1", "best_prediction_2", "best_prediction_3"], "error_bands": True, #"y_scale": (10**(-8),5*10**(-3)), }, {"data_ind": ["%s_central_scale_proVBFH_NLO_ratio"%proVBF_setup, "%s_central_scale_proVBFH_LO_ratio"%proVBF_setup, "%s_central_scale_proVBFH_NNLO_ratio"%proVBF_setup, "%s_proVBFH_NNLO_nonfact_ratio"%proVBF_setup], #"reference_ind": "bigstxs_postprocessed_proVBFH_NLO", "error_bands": True, "horizontal_line": False, #"plot_label_for_reference": True, "label": "ratio to\nVBF NLO QCD", "y_scale": (0.4,1.53), }, {"data_ind": ["full_nlo_ew_photon_only", "full_nlo_ew_w_photon"], "error_bands": False, "horizontal_line": True, "label": "Full\n NLO EW/LO", "y_scale": (0.65,1.02), }, {"data_ind": ["MoCaNLO_NLO_QCD_ratio", "MoCaNLO_LO_QCD_ratio"], "error_bands": True, "label": "Full/VBF", "horizontal_line": True, "y_scale": (0.0000,6.3000), }, {"data_ind": ["MoCaNLO_loop_induced", "MoCaNLO_loop_induced_2"],#, "MoCaNLO_ggF_LO"], "reference_ind": "MoCaNLO_LO", "horizontal_line": False, "label": "Loop-induced/LO", "y_scale": (-0.005,0.7000), "error_bands": True, }, {"data_ind": ["MoCaNLO_ggF_LO"],#, "MoCaNLO_ggF_LO"], "reference_ind": "MoCaNLO_LO", "horizontal_line": False, "label": "ggF/LO", "y_scale": (-0.0,10.5000), "error_bands": True, } ] settings = default_settings figsize = default_figsize if key_to_plot == "histogram_Mjj_pTH_truth_Mjj_vs_pTH": small_ticks = ["[300,500)","[500,700)","[700,900)","[900,1100)","[1100,$\\infty$)"] medium_ticks = ["$p_T^H \\in$[0,80)","$p_T^H \\in$[80,120)","$p_T^H \\in$[120,260)","$p_T^H \\in$[260,500)","$p_T^H \\in$[500,850)","$p_T^H \\in$[850,$\\infty$)"] settings = { "complex_indexes": { "small_ticks": small_ticks, "small_ticks_rotation": 65, "medium_ticks": medium_ticks, "medium_ticks_labelpad": 90, "big_ticks_labelpad": 120, "big_ticks": [""], #"v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "big_ticks_xlines": [], }, "x_param_labelpad": -40, "x_param_label": "$M_{jj}$", "x_param_position": 1,#-0.03, "left_title": r"pp $\longrightarrow$ H + 2j", "right_title": r"$\sqrt{S} = 13.6\,\mathrm{TeV}$", } elif key_to_plot == "histogram_Mjj_Dphijj_truth_Mjj_vs_Dphijj": reference_list[0]["label"]="$\dfrac{d^2 \sigma}{d M_{jj} d \Delta \phi_{jj}} $ [pb]" #reference_list[0]["y_scale"]=10**(-8),5*10**(-3) #QCD correction reference_list[1]["y_scale"]=(0.9411,1.4) #EW correction reference_list[2]["y_scale"]=(0.87,1.05) #Full vs VBF reference_list[3]["y_scale"]=(0.85,1.15) #loop reference_list[4]["y_scale"]=(-0.0050,0.1500) #ggF reference_list[5]["y_scale"]=(-0.0,3.500) small_ticks = ["[300,500)","[500,700)","[700,900)","[900,1100)","[1100,$\\infty$)"] medium_ticks = ["$\\Delta \\phi_{jj} \\in[0,\\pi/4)$","$\\Delta \\phi_{jj} \\in[\\pi/4,\\pi/2)$", "$\\Delta \\phi_{jj} \\in[\\pi/2,3\\pi/4)$", "$\\Delta \\phi_{jj} \\in[3\\pi/4,\\pi)$"] settings = { "complex_indexes": { "small_ticks": small_ticks, "small_ticks_rotation": 65, "medium_ticks": medium_ticks, "medium_ticks_labelpad": 90, "big_ticks_labelpad": 120, "big_ticks": [""], #"v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "big_ticks_xlines": [], }, "x_param_labelpad": -40, "x_param_label": "$M_{jj}$", "x_param_position": 1,#-0.03, "left_title": r"pp $\longrightarrow$ H + 2j", "right_title": r"$\sqrt{S} = 13.6\,\mathrm{TeV}$", } elif key_to_plot == "histogram_Mjj_Dyjj_truth_Mjj_vs_Dyjj": reference_list[0]["label"]="$\dfrac{d^2 \sigma}{d M_{jj} d \Delta y_{jj}} $ [pb]" #reference_list[0]["y_scale"]=10**(-8),0.5*10**(-2) #QCD correction reference_list[1]["y_scale"]=(0.9,1.6) #EW correction reference_list[2]["y_scale"]=(0.85,1.15) #Full vs VBF reference_list[3]["y_scale"]=(0.85,1.15) #loop reference_list[4]["y_scale"]=(-0.0050,0.4500) #ggF reference_list[5]["y_scale"]=(-0.0,4.500) small_ticks = ["[300,500)","[500,700)","[700,900)","[900,1100)","[1100,$\\infty$)"] medium_ticks =["$\\Delta y_{jj} \\in[2,4)$","$\\Delta y_{jj} \\in[4,5)$","$\\Delta y_{jj} \\in[5,6)$" ,"$\\Delta y_{jj} \\in[6,7)$","$\\Delta y_{jj} \\in[7,\\infty)$"] settings = { "complex_indexes": { "small_ticks": small_ticks, "small_ticks_rotation": 65, "medium_ticks": medium_ticks, "medium_ticks_labelpad": 90, "big_ticks_labelpad": 120, "big_ticks": [""], #"v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "big_ticks_xlines": [], }, "x_param_labelpad": -40, "x_param_label": "$M_{jj}$", "x_param_position": 1,#-0.03, "left_title": r"pp $\longrightarrow$ H + 2j", "right_title": r"$\sqrt{S} = 13.6\,\mathrm{TeV}$", } elif key_to_plot == "histogram_pTH_Dyjj_truth_pTH_vs_Dyjj": reference_list[0]["label"]="$\dfrac{d^2 \sigma}{d p_{TH} d \Delta y_{jj}} $ [pb]" #reference_list[0]["y_scale"]=10**(-60),2*10**(-2) #QCD correction reference_list[1]["y_scale"]=(0.8,1.6) #EW correction reference_list[2]["y_scale"]=(0.65,1.1) #Full vs VBF reference_list[3]["y_scale"]=(0.9,1.8) #loop reference_list[4]["y_scale"]=(-0.0050,0.300) #ggF reference_list[5]["y_scale"]=(-0.0,4.00) small_ticks = ["[0,80)","[80,120)","[120,260)","[260,500)","[500,850)","[850,$\\infty$)"] medium_ticks =["$\\Delta y_{jj} \\in[2,4)$","$\\Delta y_{jj} \\in[4,5)$","$\\Delta y_{jj} \\in[5,6)$" ,"$\\Delta y_{jj} \\in[6,7)$","$\\Delta y_{jj} \\in[7,\\infty)$"] settings = { "complex_indexes": { "small_ticks": small_ticks, "small_ticks_rotation": 65, "medium_ticks": medium_ticks, "medium_ticks_labelpad": 90, "big_ticks_labelpad": 120, "big_ticks": [""], #"v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "v_lines": [i + 1 for i in range(len(small_ticks), Nbins+1, len(small_ticks))], "big_ticks_xlines": [], }, "x_param_labelpad": -40, "x_param_label": "$p_{TH}$", "x_param_position": 1,#-0.03, "left_title": r"pp $\longrightarrow$ H + 2j", "right_title": r"$\sqrt{S} = 13.6\,\mathrm{TeV}$", } elif key_to_plot == "histogram_Mjj_pTH_STSX_truth_Mjj_vs_pTH_STSX": #reference_list[0]["y_scale"]=(10**(-2),2) reference_list[0]["label"]="$\dfrac{d^3 \sigma}{d M_{jj} d p_{TH} d \Delta \phi_{jj}} $ [pb]" #QCD correction reference_list[1]["y_scale"]=(0.92,1.5) #EW correction reference_list[2]["y_scale"]=(0.8,1.2) reference_list[2].pop("y_scale") #Full vs VBF reference_list[3]["y_scale"]=(0.9,1.6) #loop reference_list[4]["y_scale"]=(-0.0030,0.30) #ggF reference_list[5]["y_scale"]=(-0.0,4.50) settings = { "complex_indexes": { "small_ticks": ["[350,700)","[700,1000)","[1000,1500)","[1500,$\\infty$)"], "small_ticks_rotation": 65, "medium_ticks": ["$\Delta \\phi_{jj} \\in[0,\\pi / 2]$", "$\Delta \\phi_{jj} \\in[\pi / 2, \pi]$", "$\Delta \\phi_{jj} \\in[0,\\pi / 2]$", "$\Delta \\phi_{jj} \\in[\pi / 2, \pi]$"], "medium_ticks_labelpad": 90, "big_ticks_labelpad": 120, "big_ticks": ["$p_T^H \\in$[0,200)","$p_T^H \\in$[200,$\\infty$)"], "v_lines": [i + 1 for i in range(4, Nbins+1, 4)], "big_ticks_xlines": [Nbins/2 + 1], }, "x_param_labelpad": -40, "x_param_label": "$M_{jj}$", "x_param_position": 1,#-0.03, "left_title": r"pp $\longrightarrow$ H + 2j", "right_title": r"$\sqrt{S} = 13.6\,\mathrm{TeV}$", } figsize = (14,14) elif key_to_plot == "histogram_Mjj_pTH_STSX_truth_Mjj_vs_pTH_STSX_Mjj_below_350": reference_list = [ {"data_ind": ["%s_central_scale_proVBFH_LO"%proVBF_setup, "best_prediction_final"],#"best_prediction_1", "best_prediction_2", "best_prediction_3"], "error_bands": True, #"y_scale": (10**(-8),2*10**(-2)), }, {"data_ind": ["MoCaNLO_LO_QCD_ratio_Mjj_smaller_350","MoCaNLO_NLO_QCD_ratio_Mjj_smaller_350"], # "reference_ind": "mocanlo_nlo_qcd_born", "error_bands": True, "horizontal_line": False, "label": "Full\n NLO QCD/LO", "y_scale": (0.9,21), }, {"data_ind": ["full_nlo_ew_photon_only", "full_nlo_ew_w_photon"], "error_bands": False, "horizontal_line": True, "label": "Full\n NLO EW/LO", #"y_scale": (0.7,1.15), }, {"data_ind": ["MoCaNLO_loop_induced", "MoCaNLO_loop_induced_2"], "reference_ind": "MoCaNLO_LO", "horizontal_line": False, "label": "Loop Induced/LO", #"y_scale": (-0.0150,0.8000), "error_bands": True, }, {"data_ind": ["MoCaNLO_ggF_LO"],#, "MoCaNLO_ggF_LO"], "reference_ind": "MoCaNLO_LO", "horizontal_line": False, "label": "ggF/LO", "y_scale": (-0.0,80.5000), "error_bands": True, } ] reference_list[0]["label"]="$\dfrac{d^2 \sigma}{d M_{jj} d \Delta \phi_{jj}} $ [pb]" #reference_list[0]["y_scale"]=(10**(-2),2) #QCD correction #reference_list[1]["y_scale"]=(0.8,1.3) #EW correction reference_list[2]["y_scale"]=(0.925,1.051) settings = { "complex_indexes": { "small_ticks": ["[0,60)","[60,120)","[120,350)"], "small_ticks_rotation": 65, "medium_ticks": ["$\Delta \\phi_{jj} \\in[0,\\pi / 2]$", "$\Delta \\phi_{jj} \\in[\pi / 2, \pi]$"], "medium_ticks_labelpad": 90, #"v_lines": [i + 1 for i in range(4, len(data["data"]["%s_central_scale_proVBFH_NLO_ratio"%proVBF_setup])+1, 4)], "big_ticks_xlines": [Nbins/2 + 1], }, "x_param_labelpad": -40, "x_param_label": "$M_{jj}$", "x_param_position": 1,#-0.03, "left_title": r"pp $\longrightarrow$ H + 2j", "right_title": r"$\sqrt{S} = 13.6\,\mathrm{TeV}$", } figsize = (12,11) else: print("Un-recognized plot type...") raise return reference_list, settings, figsize def transform_mocanlo_STXS_data(list_mocanlo_data_raw, proVBF_setup, top_data_dict_orig_bin_width, merge_pTHjj, merge_DPhijj_pTHjj, mjj_less_350=False, denormalise_flag=False): #scale by bin width: mocanlo_data_list = [] for mocanlo_data in list_mocanlo_data_raw: mocanlo=deepcopy(mocanlo_data) #if is STXS, need to remove the first bin from MoCaNLO: if proVBF_setup == "stxs": mocanlo["data"] = mocanlo["data"][1:] # 29.01.2025 Denormalise is now in MoCaNLO reader if not mjj_less_350 and not denormalise_flag: mocanlo["data"] = [(mocanlo_entry[0], mocanlo_entry[1], mocanlo_entry[2] / float(proVBF_entry[1]-proVBF_entry[0]), mocanlo_entry[3]/float(proVBF_entry[1]-proVBF_entry[0]), mocanlo_entry[4]/float(proVBF_entry[1]-proVBF_entry[0])) for mocanlo_entry,proVBF_entry in zip(mocanlo['data'],top_data_dict_orig_bin_width["%s_central_scale_proVBFH_LO"%proVBF_setup]) ] if proVBF_setup == "stxs" and merge_pTHjj and not mjj_less_350: full_data_list = mocanlo["data"] #print(full_data_list) new_data_list = [] for i in range(16): new_bin = (float(i+1),float(i+2), full_data_list[i][2]+full_data_list[i+16][2], np.sqrt(full_data_list[i][3]**2+full_data_list[i+16][3]**2), np.sqrt(full_data_list[i][4]**2+full_data_list[i+16][4]**2)) new_data_list.append(new_bin) mocanlo["data"]=new_data_list elif proVBF_setup == "stxs" and merge_DPhijj_pTHjj and not mjj_less_350: full_data_list = mocanlo["data"] #print(full_data_list) new_data_list = [] for i in range(8): print(i, i + 16) new_bin = (float(i+1),float(i+2), full_data_list[i][2]+full_data_list[i+8][2]+full_data_list[i+16][2]+full_data_list[i+24][2], np.sqrt(full_data_list[i][3]**2+full_data_list[i+8][3]**2+full_data_list[i+16][3]**2+full_data_list[i+24][3]**2), np.sqrt(full_data_list[i][4]**2+full_data_list[i+8][4]**2+full_data_list[i+16][4]**2+full_data_list[i+24][4]**2)) new_data_list.append(new_bin) mocanlo["data"]=new_data_list elif mjj_less_350 and merge_pTHjj: full_data_list = mocanlo["data"] #print(full_data_list) new_data_list = [] for i in range(32, 35): new_bin = (float(i - 31),float(i - 30), full_data_list[i][2]+full_data_list[i+3][2], np.sqrt(full_data_list[i][3]**2+full_data_list[i+3][3]**2), np.sqrt(full_data_list[i][4]**2+full_data_list[i+3][4]**2)) new_data_list.append(new_bin) for i in range(38, 41): #print(i, i + 3) new_bin = (float(i - 34),float(i - 33), full_data_list[i][2]+full_data_list[i+3][2], np.sqrt(full_data_list[i][3]**2+full_data_list[i+3][3]**2), np.sqrt(full_data_list[i][4]**2+full_data_list[i+3][4]**2)) new_data_list.append(new_bin) mocanlo["data"]=new_data_list elif mjj_less_350 and not merge_pTHjj: full_data_list = mocanlo["data"] #print(full_data_list) new_data_list = [] for i in range(32, 44): new_bin = (float(i - 31),float(i - 30), full_data_list[i][2], full_data_list[i][3], full_data_list[i][4]) new_data_list.append(new_bin) mocanlo["data"]=new_data_list mocanlo_data_list.append(mocanlo) return mocanlo_data_list def plot_fo_2d(provbf_mapping: Dict, data_abs_path: str, save_path: str, reference_list=default_reference_list, label_mapping=default_label_mapping, color_mapping=default_color_mapping, settings=default_settings, denormalise_flag=False, separated_subplots_flag=False): # choose from available list of plots for key_to_plot_1 in keys_to_plot.keys(): keys_to_plot_2 = keys_to_plot[key_to_plot_1] Mjj_below_350 = False merge_DPhijj_pTHjj = True if key_to_plot_1 == "histogram_Mjj_pTH_STSX_truth_Mjj_vs_pTH_STSX": setup = stxs_setup proVBF_setup = "stxs" Mjj_below_350 = False elif key_to_plot_1 == "histogram_Mjj_pTH_STSX_truth_Mjj_vs_pTH_STSX_Mjj_below_350": setup = stxs_setup proVBF_setup = "stxs" Mjj_below_350 = True else: setup = fiducial_setups proVBF_setup = "histo" # do pTj 20/30 for fiducial 2D plots for setup_pTj in setup.keys(): setup2a2b = setup[setup_pTj] if setup_pTj == "stxs": keys_to_plot_PTH20 = keys_to_plot_2 else: keys_to_plot_PTH20 = [x+setup_pTj for x in keys_to_plot_2] #print(keys_to_plot_PTH20) print("producing {} under setup {}".format(key_to_plot_1, setup2a2b)) data={} data["data"] = {} top_data_dict_orig_bin_width={} reader = TopReader(mapping=provbf_mapping, denormalise_flag=denormalise_flag) if not Mjj_below_350: # Getting Mjj-ptH data from top reader provbf_histo_list0 = list(filter(lambda x: "proVBFH" in x,dir(reader))) if setup2a2b=="STXS": provbf_histo_list = list(filter(lambda x: x.startswith("stxs_"), provbf_histo_list0)) #print(provbf_histo_list) else: provbf_histo_list = list(filter(lambda x: "histo" in x, provbf_histo_list0)) top_data_dict={} #iterate through different LO/NLO/NNLO histo type for hist_type in provbf_histo_list: full_hist_2D_unrolled = [] full_hist_2D_unrolled_orig_bin_width = [] #get each of the Mjj hist in different pTH bins #print(keys_to_plot_PTH20) for itr_pTHbin, pTH_bin in enumerate(keys_to_plot_PTH20): #get the data of try: pTH_bin_data = list(filter(lambda x: pTH_bin in x["index"],getattr(reader, hist_type)))[0]["data"] except: pass for mjj_bin_count in range(0,len(pTH_bin_data)): #print(hist_type, pTH_bin_data[mjj_bin_count]) list_pTH_bin_data = list(pTH_bin_data[mjj_bin_count]) list_pTH_bin_data[0] = float(itr_pTHbin*len(pTH_bin_data)+mjj_bin_count+1) list_pTH_bin_data[1] = float(itr_pTHbin*len(pTH_bin_data)+mjj_bin_count+2) #list_pTH_bin_data.append(list_pTH_bin_data[-1]) # TODO WHY??? tuple_pTH_bin_data = tuple(list_pTH_bin_data) full_hist_2D_unrolled.append(tuple_pTH_bin_data) list_pTH_bin_data_orig_bin_width = list(pTH_bin_data[mjj_bin_count]) tuple_pTH_bin_data_orig_bin_width = tuple(list_pTH_bin_data_orig_bin_width) full_hist_2D_unrolled_orig_bin_width.append(tuple_pTH_bin_data_orig_bin_width) top_data_dict[hist_type] = full_hist_2D_unrolled top_data_dict_orig_bin_width[hist_type] = full_hist_2D_unrolled_orig_bin_width ##for STXS - merge DPhijj and pTHjj proVBF samples : if proVBF_setup == "stxs" and merge_DPhijj_pTHjj and not merge_pTHjj: import numpy as np new_data_dict = {} for proVBF_sample in top_data_dict.keys(): full_data_list = top_data_dict[proVBF_sample] new_data_dict[proVBF_sample] = [] for i in range(8): new_bin = (float(i+1),float(i+2), full_data_list[i][2]+full_data_list[i+8][2]+full_data_list[i+16][2]+full_data_list[i+24][2], np.sqrt(full_data_list[i][3]**2+full_data_list[i+8][3]**2+full_data_list[i+16][3]**2+full_data_list[i+24][3]**2), np.sqrt(full_data_list[i][4]**2+full_data_list[i+8][4]**2+full_data_list[i+16][4]**2+full_data_list[i+24][4]**2)) new_data_dict[proVBF_sample].append(new_bin) top_data_dict = new_data_dict elif merge_pTHjj and proVBF_setup == "stxs": import numpy as np new_data_dict = {} for proVBF_sample in top_data_dict.keys(): full_data_list = top_data_dict[proVBF_sample] new_data_dict[proVBF_sample] = [] for i in range(16): new_bin = (float(i+1),float(i+2), full_data_list[i][2]+full_data_list[i+16][2], np.sqrt(full_data_list[i][3]**2+full_data_list[i+16][3]**2), np.sqrt(full_data_list[i][4]**2+full_data_list[i+16][4]**2)) new_data_dict[proVBF_sample].append(new_bin) top_data_dict = new_data_dict data["data"] = top_data_dict key_to_plot_1_copy = key_to_plot_1 if key_to_plot_1 == "histogram_Mjj_pTH_STSX_truth_Mjj_vs_pTH_STSX_Mjj_below_350": key_to_plot_1_copy = "histogram_Mjj_pTH_STSX_truth_Mjj_vs_pTH_STSX" red = MocanloReader() # mocanlo_nlo_ew_no_photon_born (mocanlo_nlo_qcd_born_0, mocanlo_nlo_qcd_nlo_0, mocanlo_nlo_ew_w_photon_born_0, mocanlo_nlo_ew_w_photon_nlo_0, mocanlo_loop_induced_0, mocanlo_loop_induced_2_0, mocanlo_nlo_ew_no_photon_nlo_0, mocanlo_nlo_ew_photon_only_0, mocanlo_ggF_born_0) = red.get_mocanlo_data_from_file_path(data_abs_path, key_to_plot_1_copy, setup2a2b, denormalise_flag) list_mocanlo_data_raw = [mocanlo_nlo_qcd_born_0, mocanlo_nlo_qcd_nlo_0, mocanlo_nlo_ew_w_photon_born_0, mocanlo_nlo_ew_w_photon_nlo_0, mocanlo_loop_induced_0, mocanlo_loop_induced_2_0, mocanlo_nlo_ew_no_photon_nlo_0, mocanlo_nlo_ew_photon_only_0,mocanlo_ggF_born_0] (mocanlo_nlo_qcd_born, mocanlo_nlo_qcd_nlo, mocanlo_nlo_ew_w_photon_born, mocanlo_nlo_ew_w_photon_nlo, mocanlo_loop_induced, mocanlo_loop_induced_2, mocanlo_nlo_ew_no_photon_nlo, mocanlo_nlo_ew_photon_only, mocanlo_ggF_born) = transform_mocanlo_STXS_data(list_mocanlo_data_raw, proVBF_setup, top_data_dict_orig_bin_width, True, merge_DPhijj_pTHjj, Mjj_below_350, denormalise_flag) data["x_param"] = mocanlo_nlo_ew_w_photon_born["x_param"] print(data["x_param"]) data["other_params"] = [] #sigma_qcd_nlo mocanlo_qcd_nlo_prediction = red.sum_mocanlo_data(mocanlo_nlo_qcd_born["data"], mocanlo_nlo_qcd_nlo["data"]) one = [(i[0], i[1], 1, 0, 0) for i in mocanlo_qcd_nlo_prediction] #(1+deltaNLO EW) delta_ew = red.sum_mocanlo_data(one, mocanlo_nlo_ew_no_photon_nlo["data"]) if not Mjj_below_350: #(1+deltaNNLO QCD) delta_vbf = red.delta_mocanlo(red.sum_mocanlo_data( data["data"]["%s_proVBFH_NNLO_nonfact"%proVBF_setup], red.minus_mocanlo_data(data["data"]["%s_central_scale_proVBFH_NNLO"%proVBF_setup],data["data"]["%s_central_scale_proVBFH_NLO"%proVBF_setup])), data["data"]["%s_central_scale_proVBFH_LO"%proVBF_setup] ) one_plus_delta_vbf = red.sum_mocanlo_data(one, delta_vbf) best_2 = red.sum_mocanlo_data(red.product_top_mocanlo(mocanlo_qcd_nlo_prediction, red.product_top_mocanlo(one_plus_delta_vbf, delta_ew)), mocanlo_nlo_ew_photon_only) else: best_2 = red.sum_mocanlo_data(red.product_top_mocanlo(mocanlo_qcd_nlo_prediction, delta_ew), mocanlo_nlo_ew_photon_only) ##Delta photon (quick and dirty rescaling) #delta_ew_photon_only_1 = mocanlo_nlo_ew_photon_only["data"] #mocanlo_nlo_ew_photon_only['data'] = [ # (entry[0], entry[1], entry[2] / 200., entry[3]/200., entry[4]/200.) # for entry in mocanlo_nlo_ew_photon_only['data'] #] #################################################################################################################### full_nlo_ew_w_photon = red.sum_mocanlo_data(one, red.delta_mocanlo(mocanlo_nlo_ew_w_photon_nlo["data"], mocanlo_nlo_ew_w_photon_born["data"])) full_nlo_ew_photon_only = red.sum_mocanlo_data(one, red.delta_mocanlo(mocanlo_nlo_ew_photon_only["data"],mocanlo_nlo_ew_w_photon_born["data"])) mocanlo_correction_and_born = red.sum_mocanlo_data(mocanlo_nlo_qcd_nlo["data"], mocanlo_nlo_qcd_born["data"]) data["data"]["MoCaNLO_LO"] = mocanlo_nlo_ew_w_photon_born["data"] #delta_ew_photon_only_1 = mocanlo_nlo_ew_photon_only["data"] #mocanlo_nlo_ew_w_photon_born['data'] = [ # (entry[0], entry[1], entry[2] / 200., entry[3]/200., entry[4]/200.) # for entry in mocanlo_nlo_ew_w_photon_born['data'] #] data["data"]["MoCaNLO_loop_induced"] = mocanlo_loop_induced["data"] data["data"]["MoCaNLO_loop_induced_2"] = mocanlo_loop_induced_2["data"] data["data"]["MoCaNLO_ggF_LO"] = mocanlo_ggF_born["data"] data["data"]["best_prediction_final"] = best_2 data["data"]["full_nlo_ew_w_photon"] = full_nlo_ew_w_photon data["data"]["full_nlo_ew_photon_only"] = full_nlo_ew_photon_only if not Mjj_below_350: data["data"]["MoCaNLO_NLO_QCD_ratio"] = red.delta_mocanlo(mocanlo_correction_and_born, data["data"]["%s_central_scale_proVBFH_NLO"%proVBF_setup]) data["data"]["MoCaNLO_LO_QCD_ratio"] = red.delta_mocanlo(mocanlo_nlo_ew_w_photon_born["data"], data["data"]["%s_central_scale_proVBFH_LO"%proVBF_setup]) data["data"]["%s_central_scale_proVBFH_LO_ratio"%proVBF_setup] = red.delta_mocanlo(data["data"]["%s_central_scale_proVBFH_LO"%proVBF_setup], data["data"]["%s_central_scale_proVBFH_NLO"%proVBF_setup]) data["data"]["%s_central_scale_proVBFH_NLO_ratio"%proVBF_setup] = red.delta_mocanlo(data["data"]["%s_central_scale_proVBFH_NLO"%proVBF_setup], data["data"]["%s_central_scale_proVBFH_NLO"%proVBF_setup]) data["data"]["%s_central_scale_proVBFH_NNLO_ratio"%proVBF_setup] = red.delta_mocanlo(data["data"]["%s_central_scale_proVBFH_NNLO"%proVBF_setup], data["data"]["%s_central_scale_proVBFH_NLO"%proVBF_setup]) data["data"]["%s_proVBFH_NNLO_nonfact_ratio"%proVBF_setup] = red.sum_mocanlo_data(one, red.delta_mocanlo(data["data"]["%s_proVBFH_NNLO_nonfact"%proVBF_setup], data["data"]["%s_central_scale_proVBFH_NLO"%proVBF_setup]) ) data["data"]["MoCaNLO_NLO_QCD_ratio_Mjj_smaller_350"] = red.delta_mocanlo(mocanlo_qcd_nlo_prediction, mocanlo_nlo_qcd_born["data"]) data["data"]["MoCaNLO_LO_QCD_ratio_Mjj_smaller_350"] = red.delta_mocanlo(mocanlo_nlo_qcd_born["data"], mocanlo_nlo_qcd_born["data"]) # reference_list is just a settings json for plotting # label_mapping is a dict to rename data["data"].keys() in the plot Nbins = len(data["data"]["best_prediction_final"]) #print("Nbins: ",Nbins) #print(data["data"]["best_prediction_final"]) label_mapping = { "full_nlo_ew_photon_only": "Photon only", "full_nlo_ew_w_photon": "NLO with photon", "best_prediction_final": "Best prediction", "MoCaNLO_NLO_EW_photon_only": "Photon only", "best_prediction_1": "without photon", "best_prediction_2" : "photon only correction", "best_prediction_3" : "EW NLO with photon - EW NLO without photon", "bigstxs_postprocessed_proVBFH_LO": "LO", "bigstxs_postprocessed_proVBFH_NLO": "NLO QCD", "STXS_proVBFH_NNLO_EW": "Best prediction", "bigstxs_postprocessed_proVBFH_LO_ratio": "LO", "bigstxs_postprocessed_proVBFH_NLO_ratio": "NLO", "bigstxs_postprocessed_proVBFH_NNLO_ratio": "NNLO", "bigstxs_postprocessed_proVBFH_NNLO": "NNLO QCD", "%s_central_scale_proVBFH_LO"%proVBF_setup: "LO", "%s_central_scale_proVBFH_NLO"%proVBF_setup: "NLO", "%s_central_scale_proVBFH_NNLO"%proVBF_setup: "NNLO", "%s_proVBFH_NNLO_nonfact_ratio"%proVBF_setup: "NF", "%s_central_scale_proVBFH_LO_ratio"%proVBF_setup: "LO", "%s_central_scale_proVBFH_NLO_ratio"%proVBF_setup: "NLO", "%s_central_scale_proVBFH_NNLO_ratio"%proVBF_setup: "NNLO", "full_nlo_qcd_correction": "NLO", "mocanlo_qcd_nlo_prediction": "NLO", "mocanlo_nlo_qcd_born": "LO", "MoCaNLO_NLO_EW": "", "MoCaNLO_NLO_QCD_ratio": "NLO", "MoCaNLO_LO_QCD_ratio": "LO", "MoCaNLO_loop_induced": "Loop-induced interference", "MoCaNLO_loop_induced_2": "Loop-induced squared", "MoCaNLO_ggF_LO": "LO ggF", "MoCaNLO_NLO_QCD_ratio_Mjj_smaller_350": "NLO", "MoCaNLO_LO_QCD_ratio_Mjj_smaller_350": "LO" } color_mapping = { "full_nlo_ew_photon_only": "#E1C16E", #Brass "full_nlo_ew_w_photon": "#6F4E37", #coffee "best_prediction_final": "#50C878", #Emerald Green "MoCaNLO_NLO_EW_photon_only": "", "best_prediction_1": "", "best_prediction_2" : "", "best_prediction_3" : "", "%s_central_scale_proVBFH_LO"%proVBF_setup: "#0096FF", #bright blue "%s_central_scale_proVBFH_NLO"%proVBF_setup: "#5F9EA0", #Cadet Blue "STXS_proVBFH_NNLO_EW": "", "%s_central_scale_proVBFH_LO_ratio"%proVBF_setup: "#0096FF", "%s_central_scale_proVBFH_NLO_ratio"%proVBF_setup: "#5F9EA0", #Cadet Blue "MoCaNLO_LO_QCD_ratio_Mjj_smaller_350": "#0096FF", "MoCaNLO_NLO_QCD_ratio_Mjj_smaller_350": "#5F9EA0", #Cadet Blue "full_nlo_qcd_correction": "#5F9EA0", #Cadet Blue "mocanlo_qcd_nlo_prediction": "#5F9EA0", "mocanlo_nlo_qcd_born": "#0096FF", "%s_central_scale_proVBFH_NNLO_ratio"%proVBF_setup: "#C30000", #"#5D3FD3",#Iris "%s_proVBFH_NNLO_nonfact_ratio"%proVBF_setup: "#F5CF9F", #light orange "MoCaNLO_NLO_EW": "", "MoCaNLO_NLO_QCD_ratio": "#988558", #Dark Tan "MoCaNLO_LO_QCD_ratio": "#CD7F32", #Bronze "MoCaNLO_loop_induced": "#708090",#Slate Gray "MoCaNLO_loop_induced_2": "#8A9A5B",#Sage Green "MoCaNLO_ggF_LO":"#F5CF9F" } reference_list, settings, figsize = plotting_setup(key_to_plot_1, proVBF_setup, Nbins) settings["separated_subplots_flag"] = separated_subplots_flag plotter = HistoPlotter(reader) #plotter.get_bigstxs_complex_histogram(data, reference_list, label_mapping, 1, 31, 0, 3, settings) fig = plotter.get_bigstxs_complex_histogram_histogramlike(data, reference_list, label_mapping, 1, Nbins+1, 0, custom_color_map=color_mapping, settings = settings, figsize=figsize)#figsize = (15,12), fig.savefig("{}/{}_{}.pdf".format(save_path, key_to_plot_1,setup2a2b), bbox_inches="tight", pad_inches=0.5)