diff --git a/benchmarks/bench_cli.py b/benchmarks/bench_cli.py index 8bc1ffdf..ca518cb4 100644 --- a/benchmarks/bench_cli.py +++ b/benchmarks/bench_cli.py @@ -113,9 +113,8 @@ def benchmark_scaffold_cli(test_empty_proj): assert "Success: All the folders are correctly configured" in res.output -def benchmark_gen_sv_cli(test_files, tmp_path, test_pdf, lhapdf_path): +def benchmark_gen_sv_cli(test_files, tmp_path): runner = CliRunner() - pdf_name = "NNPDF40_nlo_as_01180" max_as = "2" nf = "5" name_grid = "ATLAS_TTB_8TEV_LJ_TRAP_norensv_fixed.pineappl.lz4" @@ -123,8 +122,6 @@ def benchmark_gen_sv_cli(test_files, tmp_path, test_pdf, lhapdf_path): new_grid_path = tmp_path / name_grid target_path = tmp_path shutil.copy(grid_path, new_grid_path) - with lhapdf_path(test_pdf): - pdf = lhapdf.mkPDF(pdf_name) res = runner.invoke( command, ["ren_sv_grid", str(new_grid_path), str(target_path), max_as, nf, "False"], diff --git a/src/pineko/kfactor.py b/src/pineko/kfactor.py index 27667ec8..9a59c2c4 100644 --- a/src/pineko/kfactor.py +++ b/src/pineko/kfactor.py @@ -9,6 +9,7 @@ from pineappl import import_only_subgrid from . import scale_variations +from . scale_variations import orders_as_tuple DEFAULT_PDF_SET = "NNPDF40_nnlo_as_01180" @@ -156,7 +157,7 @@ def construct_new_order(grid, order, order_to_update, central_kfactor, alphas): alphas """ # extract the relevant order to rescale from the grid for each lumi and bin - grid_orders = [order.as_tuple() for order in grid.orders()] + grid_orders = orders_as_tuple(grid) new_grid = scale_variations.initialize_new_grid(grid, order_to_update) orginal_order_index = grid_orders.index(order) @@ -215,14 +216,15 @@ def do_it( order_exists: bool True if the order to update is already present """ - grid_orders = [order.as_tuple() for order in grid.orders()] + grid_orders = orders_as_tuple(grid) # remove not necessary orders - # NOTE: eventual QED corrections are not supported + # NOTE: eventual QED kfactors are not supported order_mask = pineappl.grid.Order.create_mask(grid.orders(), pto_to_update, 0, True) grid_orders_filtered = list(np.array(grid_orders)[order_mask]) grid_orders_filtered.sort(key=scale_variations.qcd) min_as = grid_orders_filtered[0][0] + # TODO: this is always going to be 0, given the mask above ... min_al = grid_orders_filtered[0][1] # the actual alpha_s order to update @@ -243,7 +245,7 @@ def do_it( max_as = grid_orders_filtered[-1][0] orders_list = [(de, min_al, 0, 0) for de in range(min_as, max_as + 1)] # create an empty grid and add the rescaled order - order_to_update = (order_to_update, grid_orders_filtered[0][1], 0, 0) + order_to_update = (order_to_update, min_al, 0, 0) new_order_grid = None for i, as_order in enumerate(orders_list): order_grid = construct_new_order( diff --git a/src/pineko/scale_variations.py b/src/pineko/scale_variations.py index 67d993b9..33f9f7fe 100644 --- a/src/pineko/scale_variations.py +++ b/src/pineko/scale_variations.py @@ -30,6 +30,8 @@ def qcd(order: OrderTuple) -> int: """Extract the QCD order from an OrderTuple.""" return order[0] +def orders_as_tuple(grid: pineappl.grid.Grid) -> list[OrderTuple]: + return [order.as_tuple() for order in grid.orders()] def ren_sv_coeffs(m, max_as, logpart, which_part, nf): """Ren_sv coefficient for the requested part. @@ -114,7 +116,7 @@ def create_svonly(grid, order, new_order, scalefactor): """Create a grid containing only the renormalization scale variations at a given order for a grid.""" new_grid = initialize_new_grid(grid, new_order) # extract the relevant order to rescale from the grid for each lumi and bin - grid_orders = [order.as_tuple() for order in grid.orders()] + grid_orders = orders_as_tuple(grid) order_index = grid_orders.index(order) for lumi_index in range(len(new_grid.lumi())): for bin_index in range(grid.bins()): @@ -137,7 +139,7 @@ def create_svonly(grid, order, new_order, scalefactor): def create_grids(gridpath, max_as, nf): """Create all the necessary scale variations grids for a certain starting grid.""" grid = pineappl.grid.Grid.read(gridpath) - grid_orders = [orde.as_tuple() for orde in grid.orders()] + grid_orders = orders_as_tuple(grid) order_mask = pineappl.grid.Order.create_mask(grid.orders(), max_as, 0, True) grid_orders_filtered = list(np.array(grid_orders)[order_mask]) grid_orders_filtered.sort(key=qcd) @@ -159,47 +161,7 @@ def create_grids(gridpath, max_as, nf): ) grid_list[to_construct_order] = list_grid_order - return grid_list, nec_orders - - -def write_grids(gridpath, grid_list): - """Write the single grids.""" - base_name = gridpath.stem.split(".pineappl")[0] - final_part = ".pineappl.lz4" - grid_paths = [] - for order in grid_list: - # For each scale variation order, if more than one grid contributes, merge them all together in a single one - if len(grid_list[order]) > 1: - for grid in grid_list[order][1:]: - tmp_path = gridpath.parent / ("tmp" + final_part) - grid.write_lz4(tmp_path) - grid_list[order][0].merge_from_file(tmp_path) - tmp_path.unlink() - new_grid_path = gridpath.parent / ( - base_name + "_" + str(order[2]) + final_part - ) # order[2] is the ren_sv order - grid_paths.append(new_grid_path) - grid_list[order][0].write_lz4(new_grid_path) - return grid_paths - - -def merge_grids( - gridpath, grid_list_path, target_path=None, nec_orders=None, order_exists=False -): - """Merge the single grids in the original.""" - grid = pineappl.grid.Grid.read(gridpath) - if target_path is None: - target_path = gridpath.parent / gridpath.name - else: - target_path = target_path / gridpath.name - if order_exists: - grid = construct_and_dump_order_exists_grid( - grid, list(nec_orders)[0] if nec_orders is not None else [] - ) - for grid_path in grid_list_path: - grid.merge_from_file(grid_path) - grid_path.unlink() - grid.write_lz4(target_path) + return grid_list def construct_and_dump_order_exists_grid(ori_grid, to_construct_order): @@ -212,13 +174,11 @@ def construct_and_dump_order_exists_grid(ori_grid, to_construct_order): to_construct_order: order to delete """ - - # TODO: maybe we can make this function simpler ?? - + # TODO: can we make this function simpler ?? bin_limits = [float(bin) for bin in range(ori_grid.bins() + 1)] lumi_grid = [pineappl.lumi.LumiEntry(mylum) for mylum in ori_grid.lumi()] subgrid_params = pineappl.subgrid.SubgridParams() - ori_grid_orders = [order.as_tuple() for order in ori_grid.orders()] + ori_grid_orders = orders_as_tuple(ori_grid) new_orders = [ pineappl.grid.Order(*ord) for ord in ori_grid_orders @@ -268,9 +228,8 @@ def compute_ren_sv_grid( ): """Generate renormalization scale variation terms for the given grid, according to the max_as.""" # First let's check if the ren_sv are already there - checkres, max_as_effective = check.contains_sv( - pineappl.grid.Grid.read(grid_path), max_as, 0, check.Scale.REN - ) + grid = pineappl.grid.Grid.read(grid_path) + checkres, max_as_effective = check.contains_sv(grid, max_as, 0, check.Scale.REN) # Usual different convention with max_as if max_as_effective == max_as and (checkres is not check.AvailableAtMax.CENTRAL): if not order_exists: @@ -279,18 +238,21 @@ def compute_ren_sv_grid( return ReturnState.ORDER_EXISTS_FAILURE if max_as_effective < max_as and checkres is check.AvailableAtMax.SCVAR: return ReturnState.MISSING_CENTRAL + # Create all the necessary grids # With respect to the usual convention here max_as is max_as-1 - max_as -= 1 - # Creating all the necessary grids - grid_list, nec_orders = create_grids(grid_path, max_as, nf) - # Writing the sv grids - sv_grids_paths = write_grids(gridpath=grid_path, grid_list=grid_list) - # Merging all together - merge_grids( - gridpath=grid_path, - grid_list_path=sv_grids_paths, - target_path=target_path, - nec_orders=nec_orders, - order_exists=order_exists, - ) + grid_list = create_grids(grid_path, max_as - 1, nf) + + # merge them + for sv_order, sv_grids in grid_list.items(): + # if the new order is there, clean the old one. + if sv_order in orders_as_tuple(grid): + grid = construct_and_dump_order_exists_grid( + grid, list(grid_list)[0] + ) + for sv_grid in sv_grids: + grid.merge(sv_grid) + # save + if target_path is None: + target_path = grid_path.parent + grid.write_lz4(target_path / grid_path.name) return ReturnState.SUCCESS