Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-02-24 09:21:51

0001 import json
0002 
0003 
0004 def dumper(obj):
0005     try:
0006         return obj.toJSON()
0007     except:
0008         return obj.__dict__
0009 
0010 
0011 class index_info(object):
0012     def __init__(self, index):
0013         self.index = index
0014         self.boundaries = []
0015         self.layers_no_approach = []
0016         self.layers_with_approach = []
0017         self.approaches = []
0018         self.passiveDiscBinningR = 0
0019         self.passiveDiscBinningPhi = 0
0020         self.passiveCylinderBinningZ = 0
0021         self.passiveCylinderBinningPhi = 0
0022         self.activeBinningRorZ = 0
0023         self.activeBinningPhi = 0
0024 
0025     def __repr__(self):
0026         return repr(
0027             (
0028                 self.index,
0029                 self.boundaries,
0030                 self.layers_no_approach,
0031                 self.layers_with_approach,
0032                 self.approaches,
0033                 self.passiveDiscBinningR,
0034                 self.passiveDiscBinningPhi,
0035                 self.passiveCylinderBinningZ,
0036                 self.passiveCylinderBinningPhi,
0037                 self.activeBinningRorZ,
0038                 self.activeBinningPhi,
0039             )
0040         )
0041 
0042 
0043 def append_index_if_missing(dictionary, name, index):
0044     if name not in dictionary:
0045         dictionary[name] = index_info(index)
0046 
0047 
0048 def extract_coords(coords, is_disc):
0049     x_values = [coords[1], coords[2]]
0050     y_values = [coords[0], coords[0]]
0051     if is_disc:
0052         x_values = [coords[2], coords[2]]
0053         y_values = [coords[1], coords[0]]
0054     return x_values, y_values
0055 
0056 
0057 def dump_geo(filename, plot, output_folder, dump_steering, steering_file):
0058     f = open(filename)
0059     data = json.load(f)
0060 
0061     # First you want to read the Volume entries and save the relevant quantities
0062     # - names
0063     index_to_names = []
0064     for entry in data["Volumes"]["entries"]:
0065         index_to_names.append(entry["value"]["NAME"])
0066 
0067     # In case volume names are not found check in the volume names with dummy values
0068     if not index_to_names:
0069         for entry in data["Surfaces"]["entries"]:
0070             if "volume" in entry:
0071                 vol = entry["volume"]
0072                 if "volume" + str(vol) not in index_to_names:
0073                     index_to_names.append("volume" + str(vol))
0074 
0075     # Stores the information to easily decide on which layers you want the material to be mapped
0076     steering_map = {}
0077 
0078     # Once you have read the volumes extensions, you read the surfaces representing layers and boundaries.
0079     # All these surfaces belong to a volume, they have therefore a volume index.
0080     # You are interested only in:
0081     # - surfaces with layer index (NO approach index):
0082     #    -- even layer indices represent "active" layers, i.e. the ones corresponding to sensitive layers
0083     #    -- odd event indices represent navigation layers
0084     # - surfaces with layer index AND approach index:
0085     #    -- indicate approach layers to "active" layers:
0086     #       e.g. it can happen that for cylinders: 1 sits in front of the active layer,
0087     #                                              2 sits in the middle of the layer,
0088     #                                              3 sits behind the layer.
0089     #                                              ...
0090     # - surfaces with boundary index (no layer index in this case):
0091     #    -- they indicate boundaries between volumes. You are interested in boundaries between volumes
0092     #       containing at least an active layer.
0093 
0094     index_to_extends_layers_bounds_cylinders = [[] for _ in range(len(index_to_names))]
0095     index_to_extends_layers_bounds_discs = [[] for _ in range(len(index_to_names))]
0096 
0097     index_to_extends_layers_cylinders = [[] for _ in range(len(index_to_names))]
0098     index_to_extends_layers_discs = [[] for _ in range(len(index_to_names))]
0099 
0100     for entry in data["Surfaces"]["entries"]:
0101         if "layer" in entry:
0102             extends = []
0103             vol = entry["volume"]
0104 
0105             if "sensitive" in entry:
0106                 continue
0107 
0108             # Filling extends with the following quantities:
0109             # for cylinders:
0110             #     radius, z min, z max, approach index, layer index
0111             # for discs:
0112             #     inner radius, outer radius, z position, approach index, layer index
0113 
0114             z_shift = 0.0
0115             if entry["value"]["transform"]["translation"] != None:
0116                 z_shift = entry["value"]["transform"]["translation"][2]
0117 
0118             approach_index = -1
0119             if "approach" in entry:
0120                 approach_index = entry["approach"]
0121 
0122             if entry["value"]["type"] == "CylinderSurface":
0123                 # radius, z min, z max, approach index, layer index
0124                 extends = [
0125                     entry["value"]["bounds"]["values"][0],
0126                     z_shift - entry["value"]["bounds"]["values"][1],
0127                     z_shift + entry["value"]["bounds"]["values"][1],
0128                     approach_index,
0129                     entry["layer"],
0130                 ]
0131                 index_to_extends_layers_cylinders[vol - 1].append(extends)
0132             elif entry["value"]["type"] == "DiscSurface":
0133                 # inner radius, outer radius, z position, approach index, layer index
0134                 extends = [
0135                     entry["value"]["bounds"]["values"][0],
0136                     entry["value"]["bounds"]["values"][1],
0137                     z_shift,
0138                     approach_index,
0139                     entry["layer"],
0140                 ]
0141                 index_to_extends_layers_discs[vol - 1].append(extends)
0142             else:
0143                 print(
0144                     f"WARNING: Processing surface with unknown type '{entry["value"]["type"]}'. Only CylinderSurface and DiscSurface are considered."
0145                 )
0146 
0147         if "boundary" in entry:
0148             extends = []
0149             vol = entry["volume"]
0150 
0151             # Filling extends with the following quantities:
0152             # for cylinders:
0153             #     radius, z min, z max, boundary index
0154             # for discs:
0155             #     inner radius, outer radius, z position, boundary index
0156 
0157             z_shift = 0.0
0158             if entry["value"]["transform"]["translation"] != None:
0159                 z_shift = entry["value"]["transform"]["translation"][2]
0160 
0161             if entry["value"]["type"] == "CylinderSurface":
0162                 extends = [
0163                     entry["value"]["bounds"]["values"][0],
0164                     z_shift - entry["value"]["bounds"]["values"][1],
0165                     z_shift + entry["value"]["bounds"]["values"][1],
0166                     entry["boundary"],
0167                 ]
0168                 index_to_extends_layers_bounds_cylinders[vol - 1].append(extends)
0169             elif entry["value"]["type"] == "DiscSurface":
0170                 extends = [
0171                     entry["value"]["bounds"]["values"][0],
0172                     entry["value"]["bounds"]["values"][1],
0173                     z_shift,
0174                     entry["boundary"],
0175                 ]
0176                 index_to_extends_layers_bounds_discs[vol - 1].append(extends)
0177             else:
0178                 print(
0179                     f"WARNING: Processing surface with unknown type '{entry["value"]["type"]}'. Only CylinderSurface and DiscSurface are considered."
0180                 )
0181 
0182     # Steering the information and collect it into an output file if needed
0183     from itertools import chain
0184 
0185     interesting_volumes = []
0186     v_index = 0
0187     is_disc = False
0188     for elements in chain(
0189         index_to_extends_layers_cylinders, index_to_extends_layers_discs
0190     ):
0191         for coords in elements:
0192             if coords[3] > 0:
0193                 continue
0194             if v_index not in interesting_volumes:
0195                 interesting_volumes.append(v_index)
0196             append_index_if_missing(steering_map, index_to_names[v_index], v_index + 1)
0197             steering_map[index_to_names[v_index]].layers_no_approach.append(coords[4])
0198         v_index = v_index + 1
0199         if not is_disc and v_index == len(index_to_extends_layers_cylinders):
0200             v_index = 0
0201             is_disc = True
0202 
0203     v_index = 0
0204     is_disc = False
0205     for elements in chain(
0206         index_to_extends_layers_bounds_cylinders, index_to_extends_layers_bounds_discs
0207     ):
0208         for coords in elements:
0209             if v_index in interesting_volumes:
0210                 append_index_if_missing(
0211                     steering_map, index_to_names[v_index], v_index + 1
0212                 )
0213                 steering_map[index_to_names[v_index]].boundaries.append(coords[3])
0214         v_index = v_index + 1
0215         if not is_disc and v_index == len(index_to_extends_layers_bounds_cylinders):
0216             v_index = 0
0217             is_disc = True
0218 
0219     v_index = 0
0220     is_disc = False
0221     for elements in chain(
0222         index_to_extends_layers_cylinders, index_to_extends_layers_discs
0223     ):
0224         for coords in elements:
0225             if coords[3] == -1:
0226                 continue
0227             if (
0228                 coords[4]
0229                 not in steering_map[index_to_names[v_index]].layers_with_approach
0230             ):
0231                 steering_map[index_to_names[v_index]].layers_with_approach.append(
0232                     coords[4]
0233                 )
0234             if coords[3] not in steering_map[index_to_names[v_index]].approaches:
0235                 steering_map[index_to_names[v_index]].approaches.append(coords[3])
0236         v_index = v_index + 1
0237         if not is_disc and v_index == len(index_to_extends_layers_bounds_cylinders):
0238             v_index = 0
0239             is_disc = True
0240 
0241     if dump_steering:
0242         output_map = {"SteeringField": steering_map}
0243         with open(steering_file, "w", encoding="utf-8") as f:
0244             json.dump(output_map, f, default=dumper, ensure_ascii=False, indent=4)
0245 
0246     # Once you have all the data in hands, you make a few nice plots to show volumes and surfaces
0247     if plot:
0248         import matplotlib.pyplot as plt
0249 
0250         plt.rcParams.update({"figure.max_open_warning": 0})
0251         from matplotlib.pyplot import cm
0252         import numpy as np
0253 
0254         color = cm.rainbow(np.linspace(0, 1, len(index_to_extends_layers_cylinders)))
0255 
0256         is_in_legend = []
0257 
0258         plt.figure(figsize=(20, 10))
0259 
0260         # Plot all layers (w/o approach layers) + boundaries
0261         v_index = 0
0262         is_disc = False
0263         for elements in chain(
0264             index_to_extends_layers_cylinders, index_to_extends_layers_discs
0265         ):
0266             for coords in elements:
0267                 # skip approach layers
0268                 if coords[3] > 0:
0269                     continue
0270                 x_values, y_values = extract_coords(coords, is_disc)
0271                 if index_to_names[v_index] not in is_in_legend:
0272                     plt.plot(
0273                         x_values,
0274                         y_values,
0275                         c=color[v_index],
0276                         label="v: " + str(v_index + 1) + ", " + index_to_names[v_index],
0277                     )
0278                     is_in_legend.append(index_to_names[v_index])
0279                 else:
0280                     plt.plot(x_values, y_values, c=color[v_index])
0281             v_index = v_index + 1
0282             if not is_disc and v_index == len(index_to_extends_layers_cylinders):
0283                 v_index = 0
0284                 is_disc = True
0285 
0286         v_index = 0
0287         is_disc = False
0288         for elements in chain(
0289             index_to_extends_layers_bounds_cylinders,
0290             index_to_extends_layers_bounds_discs,
0291         ):
0292             for coords in elements:
0293                 if v_index in interesting_volumes:
0294                     x_values, y_values = extract_coords(coords, is_disc)
0295                     plt.plot(x_values, y_values, c=color[v_index])
0296             v_index = v_index + 1
0297             if not is_disc and v_index == len(index_to_extends_layers_bounds_cylinders):
0298                 v_index = 0
0299                 is_disc = True
0300 
0301         plt.xlabel("z [mm]")
0302         plt.ylabel("R [mm]")
0303         plt.title("Volumes and Layers (no approach layers)")
0304         plt.legend(loc="center left", bbox_to_anchor=(1, 0.5))
0305         plt.savefig(output_folder + "/volumes_and_layers.png")
0306 
0307         # Plot each volume: layers + approach layers
0308         v_index = 0
0309         is_disc = False
0310         approach_colors = ["black", "blue", "red", "green", "orange", "purple", "pink"]
0311         for elements in chain(
0312             index_to_extends_layers_cylinders, index_to_extends_layers_discs
0313         ):
0314             l_index = 0
0315             if not elements:
0316                 v_index = v_index + 1
0317                 continue
0318             plt.figure(figsize=(20, 10))
0319             color_layers = cm.rainbow(np.linspace(0, 1, len(elements)))
0320             has_elements = False
0321             for coords in elements:
0322                 if coords[3] > 0:
0323                     continue
0324                 has_elements = True
0325                 x_values, y_values = extract_coords(coords, is_disc)
0326                 plt.plot(
0327                     x_values,
0328                     y_values,
0329                     c=color_layers[l_index],
0330                     label="l: " + str(coords[4]),
0331                 )
0332                 l_index = l_index + 1
0333 
0334                 a_is_disc = False
0335                 count = 0
0336                 for a_coords in chain(
0337                     index_to_extends_layers_cylinders[v_index],
0338                     index_to_extends_layers_discs[v_index],
0339                 ):
0340                     if a_coords[4] == coords[4] and a_coords[3] > 0:
0341                         ax_values, ay_values = extract_coords(a_coords, a_is_disc)
0342                         plt.plot(
0343                             ax_values,
0344                             ay_values,
0345                             linestyle=(0, (5, 10)),
0346                             c=approach_colors[a_coords[3]],
0347                             label="l: " + str(coords[4]) + ", a: " + str(a_coords[3]),
0348                         )
0349                     count = count + 1
0350                     if count == len(index_to_extends_layers_cylinders[v_index]):
0351                         a_is_disc = True
0352 
0353             v_index = v_index + 1
0354             if not is_disc and v_index == len(index_to_extends_layers_cylinders):
0355                 v_index = 0
0356                 is_disc = True
0357 
0358             if not has_elements:
0359                 continue
0360 
0361             plt.xlabel("z [mm]")
0362             plt.ylabel("R [mm]")
0363             plt.title(index_to_names[v_index - 1])
0364             plt.legend(loc="center left", bbox_to_anchor=(1, 0.5))
0365             plt.savefig(output_folder + "/layers_for_volume_" + str(v_index) + ".png")
0366 
0367         plt.figure(figsize=(20, 10))
0368 
0369         # Plot boundaries
0370         v_index = 0
0371         is_disc = False
0372         for elements in chain(
0373             index_to_extends_layers_bounds_cylinders,
0374             index_to_extends_layers_bounds_discs,
0375         ):
0376             for coords in elements:
0377                 x_values, y_values = extract_coords(coords, is_disc)
0378                 if v_index in interesting_volumes:
0379                     plt.plot(
0380                         x_values,
0381                         y_values,
0382                         linestyle="--",
0383                         c=color[v_index],
0384                         label="v: " + str(v_index + 1) + ", b: " + str(coords[3]),
0385                     )
0386             v_index = v_index + 1
0387             if not is_disc and v_index == len(index_to_extends_layers_bounds_cylinders):
0388                 v_index = 0
0389                 is_disc = True
0390 
0391         plt.xlabel("z [mm]")
0392         plt.ylabel("R [mm]")
0393         plt.title("Boundary surfaces")
0394         plt.legend(loc="center left", bbox_to_anchor=(1, 0.5))
0395         plt.savefig(output_folder + "/boundaries.png")
0396 
0397         plt.figure(figsize=(20, 10))
0398 
0399         # Plot all approach layers
0400         v_index = 0
0401         is_disc = False
0402         add_to_legend = []
0403         for elements in chain(
0404             index_to_extends_layers_cylinders, index_to_extends_layers_discs
0405         ):
0406             if not elements:
0407                 v_index = v_index + 1
0408                 continue
0409             for coords in elements:
0410                 if coords[3] == -1:
0411                     continue
0412                 x_values, y_values = extract_coords(coords, is_disc)
0413                 if coords[3] not in add_to_legend:
0414                     plt.plot(
0415                         x_values,
0416                         y_values,
0417                         c=approach_colors[coords[3]],
0418                         linestyle="--",
0419                         label="approach index = " + str(coords[3]),
0420                     )
0421                     add_to_legend.append(coords[3])
0422                 else:
0423                     plt.plot(
0424                         x_values, y_values, c=approach_colors[coords[3]], linestyle="--"
0425                     )
0426             v_index = v_index + 1
0427             if not is_disc and v_index == len(index_to_extends_layers_bounds_cylinders):
0428                 v_index = 0
0429                 is_disc = True
0430 
0431         plt.xlabel("z [mm]")
0432         plt.ylabel("R [mm]")
0433         plt.title("Approach layers")
0434         plt.legend(loc="center left", bbox_to_anchor=(1, 0.5))
0435         plt.savefig(output_folder + "/approach_layers.png")
0436 
0437 
0438 def read_and_modify(filename, plot, output_folder, steering_file, output_file):
0439     with open(filename) as f:
0440         data = json.load(f)
0441 
0442     with open(steering_file) as f:
0443         full_data = json.load(f)
0444 
0445     layer_data = full_data["SteeringField"]
0446 
0447     index_to_names = []
0448 
0449     # Allows to dump the binning defined for the material layers
0450     dump_binning_for_material = False
0451 
0452     # Allows to check which layers are configured to carry material
0453     check_material_layers = True
0454 
0455     # First you want to read the Volume entries and save the relevant quantities
0456     # - names
0457     for entry in data["Volumes"]["entries"]:
0458         index_to_names.append(entry["value"]["NAME"])
0459 
0460     # In case volume names are not found check in the volume names with dummy values
0461     if not index_to_names:
0462         for entry in data["Surfaces"]["entries"]:
0463             if "volume" in entry:
0464                 vol = entry["volume"]
0465                 if "volume" + str(vol) not in index_to_names:
0466                     index_to_names.append("volume" + str(vol))
0467 
0468     for entry in data["Surfaces"]["entries"]:
0469         if "layer" in entry:
0470             vol = entry["volume"]
0471 
0472             if index_to_names[vol - 1] in layer_data:
0473                 if "approach" in entry:
0474                     if (
0475                         entry["layer"]
0476                         in layer_data[index_to_names[vol - 1]]["layers_with_approach"]
0477                         and entry["approach"]
0478                         in layer_data[index_to_names[vol - 1]]["approaches"]
0479                     ):
0480                         entry["value"]["material"]["mapMaterial"] = True
0481                 else:
0482                     if (
0483                         entry["layer"]
0484                         in layer_data[index_to_names[vol - 1]]["layers_no_approach"]
0485                     ):
0486                         entry["value"]["material"]["mapMaterial"] = True
0487 
0488                 if entry["value"]["material"]["mapMaterial"]:
0489                     for val in entry["value"]["material"]["binUtility"]["binningdata"]:
0490                         if val["value"] == "binZ" or val["value"] == "binR":
0491                             val["bins"] = layer_data[index_to_names[vol - 1]][
0492                                 "activeBinningRorZ"
0493                             ]
0494                         else:
0495                             val["bins"] = layer_data[index_to_names[vol - 1]][
0496                                 "activeBinningPhi"
0497                             ]
0498                         if val["bins"] == 0:
0499                             print(
0500                                 "ERROR!!! Using binning value == 0! Check you input for",
0501                                 index_to_names[vol - 1],
0502                             )
0503                             return
0504 
0505             approach_index = "None"
0506             if "approach" in entry:
0507                 approach_index = entry["approach"]
0508 
0509             if dump_binning_for_material and entry["value"]["material"]["mapMaterial"]:
0510                 print(
0511                     "Volume: ",
0512                     entry["volume"],
0513                     index_to_names[vol - 1],
0514                     " - Layer: ",
0515                     entry["layer"],
0516                     " - Approach:",
0517                     approach_index,
0518                 )
0519                 for val in entry["value"]["material"]["binUtility"]["binningdata"]:
0520                     print("-->", val["value"], ": ", val["bins"])
0521 
0522         if "boundary" in entry:
0523             extends = []
0524             vol = entry["volume"]
0525 
0526             if (
0527                 index_to_names[vol - 1] in layer_data
0528                 and entry["boundary"]
0529                 in layer_data[index_to_names[vol - 1]]["boundaries"]
0530             ):
0531                 entry["value"]["material"]["mapMaterial"] = True
0532                 for val in entry["value"]["material"]["binUtility"]["binningdata"]:
0533                     if entry["value"]["type"] == "CylinderSurface":
0534                         if val["value"] == "binZ":
0535                             val["bins"] = layer_data[index_to_names[vol - 1]][
0536                                 "passiveCylinderBinningZ"
0537                             ]
0538                         else:
0539                             val["bins"] = layer_data[index_to_names[vol - 1]][
0540                                 "passiveCylinderBinningPhi"
0541                             ]
0542                     elif entry["value"]["type"] == "DiscSurface":
0543                         if val["value"] == "binR":
0544                             val["bins"] = layer_data[index_to_names[vol - 1]][
0545                                 "passiveDiscBinningR"
0546                             ]
0547                         else:
0548                             val["bins"] = layer_data[index_to_names[vol - 1]][
0549                                 "passiveDiscBinningPhi"
0550                             ]
0551                     else:
0552                         print(
0553                             f"WARNING: Processing surface with unknown type '{entry["value"]["type"]}. Only CylinderSurface and DiscSurface are considered."
0554                         )
0555                     if val["bins"] == 0:
0556                         print(
0557                             "ERROR!!! Using binning value == 0! Check you input for",
0558                             index_to_names[vol - 1],
0559                         )
0560                         return
0561 
0562             if dump_binning_for_material and entry["value"]["material"]["mapMaterial"]:
0563                 print(
0564                     "Volume: ",
0565                     entry["volume"],
0566                     index_to_names[vol - 1],
0567                     " - Boundary:",
0568                     entry["boundary"],
0569                 )
0570                 for val in entry["value"]["material"]["binUtility"]["binningdata"]:
0571                     print("-->", val["value"], ": ", val["bins"])
0572 
0573     # Once you have all the data in hands, you make a few nice plots to show volumes and surfaces
0574 
0575     with open(output_file, "w", encoding="utf-8") as f:
0576         json.dump(data, f, ensure_ascii=False, indent=4)
0577 
0578     if plot and check_material_layers:
0579         import matplotlib.pyplot as plt
0580         from matplotlib.pyplot import cm
0581         import numpy as np
0582 
0583         plt.figure(figsize=(20, 10))
0584 
0585         material_layer_cylinders = [[] for _ in range(len(index_to_names))]
0586         material_layer_discs = [[] for _ in range(len(index_to_names))]
0587 
0588         material_approach_cylinders = [[] for _ in range(len(index_to_names))]
0589         material_approach_discs = [[] for _ in range(len(index_to_names))]
0590 
0591         material_boundary_cylinders = [[] for _ in range(len(index_to_names))]
0592         material_boundary_discs = [[] for _ in range(len(index_to_names))]
0593 
0594         for entry in data["Surfaces"]["entries"]:
0595             if not entry["value"]["material"]["mapMaterial"]:
0596                 continue
0597 
0598             z_shift = 0.0
0599             if entry["value"]["transform"]["translation"] != None:
0600                 z_shift = entry["value"]["transform"]["translation"][2]
0601 
0602             if "layer" in entry:
0603                 extends = []
0604                 vol = entry["volume"]
0605 
0606                 if entry["value"]["type"] == "CylinderSurface":
0607                     extends = [
0608                         entry["value"]["bounds"]["values"][0],
0609                         z_shift - entry["value"]["bounds"]["values"][1],
0610                         z_shift + entry["value"]["bounds"]["values"][1],
0611                     ]
0612                     if "approach" in entry:
0613                         material_approach_cylinders[vol - 1].append(extends)
0614                     else:
0615                         material_layer_cylinders[vol - 1].append(extends)
0616 
0617                 elif entry["value"]["type"] == "DiscSurface":
0618                     extends = [
0619                         entry["value"]["bounds"]["values"][0],
0620                         entry["value"]["bounds"]["values"][1],
0621                         z_shift,
0622                     ]
0623                     if "approach" in entry:
0624                         material_approach_discs[vol - 1].append(extends)
0625                     else:
0626                         material_layer_discs[vol - 1].append(extends)
0627                 else:
0628                     print(
0629                         f"WARNING: Processing surface with unknown type '{entry["value"]["type"]}'. Only CylinderSurface and DiscSurface are considered."
0630                     )
0631 
0632             if "boundary" in entry:
0633                 extends = []
0634                 vol = entry["volume"]
0635 
0636                 if entry["value"]["type"] == "CylinderSurface":
0637                     extends = [
0638                         entry["value"]["bounds"]["values"][0],
0639                         z_shift - entry["value"]["bounds"]["values"][1],
0640                         z_shift + entry["value"]["bounds"]["values"][1],
0641                     ]
0642                     material_boundary_cylinders[vol - 1].append(extends)
0643 
0644                 elif entry["value"]["type"] == "DiscSurface":
0645                     extends = [
0646                         entry["value"]["bounds"]["values"][0],
0647                         entry["value"]["bounds"]["values"][1],
0648                         z_shift,
0649                     ]
0650                     material_boundary_discs[vol - 1].append(extends)
0651                 else:
0652                     print(
0653                         "WARNING: Processing surface with unknown type. Only CylinderSurface and DiscSurface are considered."
0654                     )
0655 
0656         from itertools import chain
0657 
0658         v_index = 0
0659         is_first = True
0660         is_disc = False
0661         for elements in chain(material_layer_cylinders, material_layer_discs):
0662             l_index = 0
0663             for coords in elements:
0664                 x_values, y_values = extract_coords(coords, is_disc)
0665                 if is_first:
0666                     plt.plot(x_values, y_values, c="black", label="layer")
0667                     is_first = False
0668                 else:
0669                     plt.plot(x_values, y_values, c="black")
0670                 l_index = l_index + 1
0671             v_index = v_index + 1
0672             if not is_disc and v_index == len(material_layer_cylinders):
0673                 is_disc = True
0674                 v_index = 0
0675 
0676         v_index = 0
0677         is_first = True
0678         is_disc = False
0679         for elements in chain(material_approach_cylinders, material_approach_discs):
0680             l_index = 0
0681             for coords in elements:
0682                 x_values, y_values = extract_coords(coords, is_disc)
0683                 if is_first:
0684                     plt.plot(x_values, y_values, c="red", label="approach")
0685                     is_first = False
0686                 else:
0687                     plt.plot(x_values, y_values, c="red")
0688                 l_index = l_index + 1
0689             v_index = v_index + 1
0690             if not is_disc and v_index == len(material_approach_cylinders):
0691                 is_disc = True
0692                 v_index = 0
0693 
0694         v_index = 0
0695         is_first = True
0696         is_disc = False
0697         for elements in chain(material_boundary_cylinders, material_boundary_discs):
0698             l_index = 0
0699             for coords in elements:
0700                 x_values, y_values = extract_coords(coords, is_disc)
0701                 if is_first:
0702                     plt.plot(x_values, y_values, c="blue", label="boundary")
0703                     is_first = False
0704                 else:
0705                     plt.plot(x_values, y_values, c="blue")
0706                 l_index = l_index + 1
0707             v_index = v_index + 1
0708             if not is_disc and v_index == len(material_boundary_cylinders):
0709                 is_disc = True
0710                 v_index = 0
0711 
0712         plt.xlabel("z [mm]")
0713         plt.ylabel("R [mm]")
0714         plt.title("Layers with material")
0715         plt.legend(loc="center left", bbox_to_anchor=(1, 0.5))
0716         plt.savefig(output_folder + "/material_layers.png")
0717 
0718 
0719 import argparse
0720 import os
0721 
0722 # Initialize parser
0723 parser = argparse.ArgumentParser()
0724 parser.add_argument(
0725     "--geometry",
0726     default="",
0727     type=str,
0728     help="Specify the input file for geometry visualisation",
0729 )
0730 parser.add_argument(
0731     "--plot",
0732     default=True,
0733     action="store_true",
0734     help="Enable plot creation for geometry visualisation (Default : True)",
0735 )
0736 parser.add_argument(
0737     "--output_folder",
0738     default="plots",
0739     type=str,
0740     help="Specify the output folder for plots (Default : plots)",
0741 )
0742 parser.add_argument(
0743     "--dump_steering",
0744     default=False,
0745     action="store_true",
0746     help="Enable production of steering file for material mapping (Default : False)",
0747 )
0748 parser.add_argument(
0749     "--edit",
0750     default=False,
0751     action="store_true",
0752     help="Enable editing of input file for creation of json for material mapping (Default : False)",
0753 )
0754 parser.add_argument(
0755     "--steering_file",
0756     default="",
0757     type=str,
0758     help="Specify the steering file to guide editing of json for material mapping",
0759 )
0760 parser.add_argument(
0761     "--output_map",
0762     default="",
0763     type=str,
0764     help="Specify the output json for material mapping",
0765 )
0766 args = parser.parse_args()
0767 
0768 print(" --- Geometry visualisation and material map fie creation --- ")
0769 print()
0770 
0771 if not args.geometry:
0772     print("Error: Missing input geometry file. Please specify --geometry.")
0773     exit()
0774 
0775 if not os.path.exists(args.geometry):
0776     print("Error: Invalid file path/name in --geometry. Please check your input!")
0777     exit()
0778 
0779 if not args.geometry.endswith(".json"):
0780     print("Error: Invalid file format in --geometry. Please check your input!")
0781     exit()
0782 
0783 print("\t parsing file : ", args.geometry)
0784 if args.plot:
0785     print("\t job configured to produce plots in output folder: ", args.output_folder)
0786     if not os.path.exists(args.output_folder):
0787         os.mkdir(args.output_folder)
0788 
0789 if args.dump_steering and args.edit:
0790     print(
0791         "Error: Wrong job configuration. --dump_steering and --edit can't be \
0792         both true at the same time."
0793     )
0794     print(
0795         "\t Decide if you want to dump the steering file OR to read an existing file for editing the geometry file."
0796     )
0797     exit()
0798 
0799 if args.dump_steering:
0800     if not args.steering_file:
0801         print("Error: Missing output steering file. Please specify --steering_file.")
0802         exit()
0803     if not args.steering_file.endswith(".json"):
0804         print("Error: Invalid file format in --steering_file. It must end with .json!")
0805         exit()
0806     print(
0807         "\t job configured to produce steering file for material mapping with name: ",
0808         args.steering_file,
0809     )
0810 
0811 if args.edit:
0812     print("\t job configured to edit the input geometry file following a steering file")
0813     if not args.steering_file:
0814         print("Error: Missing input steering file. Please specify --steering_file.")
0815         exit()
0816     if not os.path.exists(args.steering_file):
0817         print(
0818             "Error: Invalid file path/name in --steering_file. Please check your input!"
0819         )
0820         exit()
0821     if not args.steering_file.endswith(".json"):
0822         print("Error: Invalid file format in --steering_file. Please check your input!")
0823         exit()
0824     if not args.output_map:
0825         print("Error: Missing output map file. Please specify --output_map.")
0826         exit()
0827     if not args.output_map.endswith(".json"):
0828         print("Error: Invalid file format in --output_map. Please check your input!")
0829         exit()
0830 
0831 
0832 if args.plot or args.dump_steering:
0833     dump_geo(
0834         args.geometry,
0835         args.plot,
0836         args.output_folder,
0837         args.dump_steering,
0838         args.steering_file,
0839     )
0840 
0841 if args.edit:
0842     read_and_modify(
0843         args.geometry,
0844         args.plot,
0845         args.output_folder,
0846         args.steering_file,
0847         args.output_map,
0848     )