Source code for smefit

# -*- coding: utf-8 -*-
[docs]class RUNNER: """ Class containing all the possible smefit functions """ def __init__(self, path): """ Init the root path of the package where tables, results plot_config and reports will be stored Parameters ---------- path : str root path """ self.root_path = path # TODO: add here a nice banner
[docs] def setup_config(self, filename): """Read yaml card and update the configuration. Parameters ---------- filename : str fit card name Returns ------- config: dict configuration dict """ import subprocess import yaml from shutil import copyfile config = {} with open(f"{self.root_path}/run_cards/{filename}.yaml") as f: config = yaml.safe_load(f) # Root path config["root_path"] = f"{self.root_path}" # Set paths to theory tables config["corrections_path"] = f"{config['root_path']}/tables/operator_res/{config['order']}/" config["theory_path"] = f"{config['root_path']}/tables/theory/" # Construct results folder subprocess.call(f"mkdir -p {config['root_path']}/results", shell=True) subprocess.call(f"mkdir -p {config['root_path']}/results/{config['resultID']}", shell=True) config["results_path"] = (f"{config['root_path']}/results/{config['resultID']}") # Copy yaml file to results folder copyfile(f"{self.root_path}/run_cards/{filename}.yaml", f"{config['results_path']}/{config['resultID']}.yaml" ) return config
[docs] def ns(self, input_card): """ Run a fit with NS given the fit name Parameters ---------- input_card : dict or str fit configuation dict or card name """ from .optimize import optimize_ns as ns from mpi4py import MPI comm = MPI.COMM_WORLD rank = comm.Get_rank() if rank == 0: if isinstance(input_card, dict): config = input_card else: config = self.setup_config(input_card) # Replica number not important for NS config.update({"rep": 1}) else: config = None config = comm.bcast(config, root=0) # Run optimizer optimize = ns.OPTIMIZE(config) optimize.run_sampling()
[docs] def individual(self, input_card): """ Run a sequence of indifidual NS fits for all the oparators set to True in the runcard Parameters ---------- input_card : str fit card name """ import copy # load general configuration base_config = self.setup_config(input_card) # loop on all the oparators for fit_op, vals in base_config["coefficients"].items(): # skip all the fixed oparators if vals["fixed"] is not False: continue individual_config = copy.deepcopy(base_config) # fix all the other oparators for op in individual_config["coefficients"]: if op != fit_op: individual_config["coefficients"][op]["fixed"] = True individual_config["coefficients"][op]["value"] = 0.0 print("Fitting the operator: ", fit_op) # run the fit self.ns(individual_config)
[docs] def mcfit(self, fit_card, scan_card=None, scan_only=False, nrep=100, first_rep=0): """ Run a fit with McFit given the input cards and the numeber of replicas. As first step it performs a `:math:\chi^2` scan then perform the MC fit of each replica. Parameters ---------- fit_card : str fit card name or scan card name for scan only scan_card : str scan card name, optional for fitting scan_only: bool perform only the chi2 scan nrep : int, number or replicas, optional for fitting first_rep: int, optional starting replica number, to be used in case you are running MCfit in parellel on a cluster """ from .optimize import chi2_profiles as chi2 from .optimize import optimize as mcfit if scan_card is None: scan_card = fit_card # chi2 scan chi2_config = self.setup_config(scan_card) chi2_config.update({"rep": 0, "crossval": False, "bootstrap": False}) chi2_prof = chi2.CHI2_PROFILES(chi2_config) chi2_prof.get_chi2_profiles() if scan_only: return if nrep <= 1: raise RuntimeError(f"nrep must be greater than 1: {nrep} given") if first_rep != 0: print(f"First replica is: {first_rep}") # MC fit fit_config = self.setup_config(fit_card) for i in range(first_rep, first_rep + nrep): # Replica fit_config.update({"rep": i}) # Run optimizer optimize = mcfit.OPTIMIZE(fit_config) optimize.run_fit()
[docs] def analysis(self, outputname, fits): """ Run the analysis given the fit ids Parameters ---------- outputname : str report directory name fits : list list of fits id """ import subprocess import numpy as np from PyPDF2 import PdfFileReader, PdfFileWriter from .analyze import analyze report_path = f"{self.root_path}/reports" dir_path = f"{report_path}/{outputname}" # Clean output folder if exists try: subprocess.call(f"rm -rf {dir_path}", shell=True) except FileNotFoundError: pass subprocess.call(f"mkdir -p {report_path}", shell=True) subprocess.call(f"mkdir -p {dir_path}", shell=True) # Initialize ANALYZE class report = analyze.ANALYZE(self.root_path, fits, outputname) # Things to include in report (single/two-parameter fits treated separately) if np.any(["SNS_" in k for k in fits]) and np.all( ["SNS_c" not in k for k in fits] ): report.write_summary() report.plot_coefficients() elif np.any(["NS_c" in k for k in fits]): report.write_summary() report.plot_ellipse() else: report.write_summary() report.write_chi2_table() report.plot_chi2() #report.write_Fisher_table() report.write_PCA_table() report.plot_coefficients() report.plot_correlations() report.plot_data_vs_theory() # Combine PDF files together into raw pdf report subprocess.call( f"gs -q -dNOPAUSE -dBATCH -sDEVICE=pdfwrite -sOutputFile={dir_path}/report_{outputname}_raw.pdf `ls -rt {dir_path}/*.pdf`", shell=True ) subprocess.call(f"mv {dir_path}/*.* {dir_path}/meta/.", shell=True) subprocess.call(f"mv {dir_path}/meta/report_*.pdf {dir_path}/.", shell=True) # Rotate PDF pages if necessary and create final report pdf_in = open(f"{dir_path}/report_{outputname}_raw.pdf", "rb") pdf_reader = PdfFileReader(pdf_in) pdf_writer = PdfFileWriter() for pagenum in range(pdf_reader.numPages): pdfpage = pdf_reader.getPage(pagenum) orientation = pdfpage.get("/Rotate") if orientation == 90: pdfpage.rotateCounterClockwise(90) pdf_writer.addPage(pdfpage) pdf_out = open(f"{dir_path}/report_{outputname}.pdf", "wb") pdf_writer.write(pdf_out) pdf_out.close() pdf_in.close() # Remove old (raw) PDF file subprocess.call(f"rm {dir_path}/report_{outputname}_raw.pdf", shell=True)