You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

110 lines
3.4 KiB

7 months ago
7 months ago
  1. # flake8: noqa
  2. from parameterisations.utils.parse_regression_coef_to_array import (
  3. parse_regression_coef_to_array,
  4. )
  5. from parameterisations.utils.fit_linear_regression_model import (
  6. fit_linear_regression_model,
  7. )
  8. import uproot
  9. import argparse
  10. from pathlib import Path
  11. def parameterise_magnet_kink(
  12. input_file: str = "data/tracking_losses_ntuple_B_default_selected.root",
  13. tree_name: str = "Selected",
  14. ) -> Path:
  15. """Function that calculates parameters for estimating the magnet kink z position.
  16. Args:
  17. input_file (str, optional): Defaults to "data/param_data_selected.root".
  18. tree_name (str, optional): Defaults to "Selected".
  19. per_layer (bool, optional): If true also calculates parameters per SciFi layer. Defaults to False.
  20. Returns:
  21. Path: Path to cpp code file.
  22. """
  23. input_tree = uproot.open({input_file: tree_name})
  24. # this is an event list of dictionaries containing awkward arrays
  25. array = input_tree.arrays()
  26. # array = allarray[
  27. # (allarray.ideal_state_770_x == allarray.ideal_state_770_x)
  28. # & (allarray.ideal_state_9410_x == allarray.ideal_state_9410_x)
  29. # & (allarray.ideal_state_10000_x == allarray.ideal_state_10000_x)
  30. # ]
  31. array["dSlope_xEndT"] = array["ideal_state_9410_tx"] - array["ideal_state_770_tx"]
  32. array["dSlope_xEndT_abs"] = abs(array["dSlope_xEndT"])
  33. array["x_EndT_abs"] = abs(array["ideal_state_9410_x"])
  34. # the magnet kink position is the point of intersection of the track tagents
  35. array["z_mag_xEndT"] = (
  36. array["ideal_state_770_x"]
  37. - array["ideal_state_9410_x"]
  38. - array["ideal_state_770_tx"] * array["ideal_state_770_z"]
  39. + array["ideal_state_9410_tx"] * array["ideal_state_9410_z"]
  40. ) / array["dSlope_xEndT"]
  41. # array = array[(array["z_mag_xEndT"] < 5500) & (array["z_mag_xEndT"] > 5000)]
  42. model_endt, poly_features_endt = fit_linear_regression_model(
  43. array,
  44. target_feat="z_mag_xEndT",
  45. features=[
  46. "ideal_state_770_tx",
  47. "dSlope_xEndT",
  48. "dSlope_xEndT_abs",
  49. "x_EndT_abs",
  50. ],
  51. keep=[
  52. "ideal_state_770_tx^2",
  53. "dSlope_xEndT^2",
  54. "dSlope_xEndT_abs",
  55. "x_EndT_abs",
  56. ],
  57. degree=2,
  58. fit_intercept=True,
  59. )
  60. cpp_ref = parse_regression_coef_to_array(
  61. model_endt,
  62. poly_features_endt,
  63. "zMagnetParamsEndT",
  64. )
  65. outpath = Path("parameterisations/result/z_mag_kink_params_electron.hpp")
  66. outpath.parent.mkdir(parents=True, exist_ok=True)
  67. with open(outpath, "w") as result:
  68. result.writelines(cpp_ref)
  69. return outpath
  70. if __name__ == "__main__":
  71. parser = argparse.ArgumentParser()
  72. parser.add_argument(
  73. "--input-file",
  74. type=str,
  75. help="Path to the input file",
  76. required=False,
  77. )
  78. parser.add_argument(
  79. "--tree-name",
  80. type=str,
  81. help="Path to the input file",
  82. required=False,
  83. )
  84. args = parser.parse_args()
  85. args_dict = {arg: val for arg, val in vars(args).items() if val is not None}
  86. outfile = parameterise_magnet_kink(**args_dict)
  87. try:
  88. import subprocess
  89. # run clang-format for nicer looking result
  90. subprocess.run(
  91. [
  92. "clang-format",
  93. "-i",
  94. f"{outfile}",
  95. ],
  96. check=True,
  97. )
  98. except:
  99. pass