Merge branch 'master' into arc_corr_and_glob_check_fix

CHANGELOG.md

@@ -14,6 +14,14 @@
   - Cleaning: Filter BPMs with NaNs
   - Cleaning: Log bad BPMs with reasons before raising errors.
 
+#### 2025-01-06 - v0.20.4 - _fsoubelet_
+
+- Fixed:
+  - Solved an issue in datetime operations occuring on `Python 3.13`.
+
+- Changed:
+  - Dropped support for `Python 3.9`.
+
 #### 2024-11-21 - v0.20.3 - _jdilly_
 
 - Fixed:

omc3/__init__.py

@@ -11,7 +11,7 @@
 __title__ = "omc3"
 __description__ = "An accelerator physics tools package for the OMC team at CERN."
 __url__ = "https://github.com/pylhc/omc3"
-__version__ = "0.20.3"
+__version__ = "0.20.4"
 __author__ = "pylhc"
 __author_email__ = "pylhc@github.com"
 __license__ = "MIT"

omc3/utils/time_tools.py

@@ -104,11 +104,12 @@ def __new__(cls, *args, **kwargs):
             dt = datetime.__new__(cls, *args, **kwargs)
 
         if dt.tzinfo is None:
-            dt = dt.replace(tzinfo=tz.tzutc())
+            if 'tzinfo' not in kwargs and len(args) < 8:  # allows forcing tz to `None`
+                dt = dt.replace(tzinfo=tz.tzutc())
 
         return datetime.__new__(cls, dt.year, dt.month, dt.day,
                                 dt.hour, dt.minute, dt.second, dt.microsecond,
-                                tzinfo=dt.tzinfo)
+                                tzinfo=dt.tzinfo, fold=dt.fold)
 
     @property
     def datetime(self):
@@ -197,7 +198,7 @@ def from_utc(cls, dt):
     @classmethod
     def from_timestamp(cls, ts):
         """Create `AccDatetime` object from timestamp."""
-        return cls(datetime.utcfromtimestamp(ts))
+        return cls(datetime.fromtimestamp(ts, tz=tz.UTC))
 
     @classmethod
     def now(cls):

pyproject.toml

@@ -24,18 +24,18 @@ authors = [
 ]
 license = "MIT"
 dynamic = ["version"]
-requires-python = ">=3.9"
+requires-python = ">=3.10"
 
 classifiers = [
     "Intended Audience :: Science/Research",
     "License :: OSI Approved :: MIT License",
     "Natural Language :: English",
     "Operating System :: OS Independent",
     "Programming Language :: Python :: 3 :: Only",
-    "Programming Language :: Python :: 3.9",
     "Programming Language :: Python :: 3.10",
     "Programming Language :: Python :: 3.11",
     "Programming Language :: Python :: 3.12",
+    "Programming Language :: Python :: 3.13",
     "Programming Language :: Python :: Implementation :: CPython",
     "Topic :: Scientific/Engineering :: Physics",
     "Topic :: Scientific/Engineering :: Visualization",
@@ -44,23 +44,20 @@ classifiers = [
     "Typing :: Typed",
 ]
 
-# TODO: when we drop python 3.9, require pytables 3.10.1 minimum and stop caring about numpy 2 compability
 dependencies = [
-  "numpy >= 1.24, < 2.0; python_version < '3.10'",  # first pytables compatible with numpy 2 is 3.10 but does not support python 3.9
-  "numpy >= 1.24; python_version >= '3.10'",  # otherwise we can use numpy 2 as on python 3.10 there is a pytables which is ok with it
+  "numpy >= 1.24",
   "scipy >= 1.10",
   "pandas >= 2.1",
   "tfs-pandas >= 3.8",
   "matplotlib >= 3.8",
-  "Pillow >= 6.2.2",  # not our dependency but older versions crash with mpl
   "generic-parser >= 1.1",
   "sdds >= 0.4",
   "optics-functions >= 0.1",
   "turn_by_turn >= 0.6",
   "uncertainties >= 3.1",
   "scikit-learn >= 1.0",
   "h5py >= 2.9",
-  "tables >= 3.9",  # TODO: require 3.10.1 minimum when it's out and we drop python 3.9 support
+  "tables >= 3.10.1",
   "requests >= 2.27",
 ]
 

tests/accuracy/test_lhc_rdts.py

@@ -0,0 +1,160 @@
+import shutil
+from pathlib import Path
+
+import numpy as np
+import pytest
+import tfs
+
+from omc3.definitions.constants import PI2
+from tests.utils.compression import decompress_model
+from tests.utils.lhc_rdts.constants import (
+    DATA_DIR,
+    MODEL_ANALYTICAL_PREFIX,
+    MODEL_NG_PREFIX,
+)
+from tests.utils.lhc_rdts.functions import (
+    get_file_suffix,
+    get_model_dir,
+    get_rdt_names,
+    get_rdt_type,
+    get_rdts_from_optics_analysis,
+    run_harpy,
+)
+
+INPUTS = Path(__file__).parent.parent / "inputs"
+
+
+@pytest.fixture(scope="module")
+def rdts_from_optics_analysis(tmp_path_factory: pytest.TempPathFactory) -> dict:
+    """Run harpy and then retrieve all the rdts from the optics analysis
+
+    This fixture decompresses the model files for both beams, runs harpy to calculate the
+    optics, and then retrieves the RDTs from the optics analysis.
+
+    Args:
+        initialise_test_paths (dict): A dictionary mapping beam number to the temporary path.
+
+    Returns:
+        dict: A dictionary mapping beam number to the RDTs calculated by OMC3.
+    """
+    dfs = {}
+    for beam in (1, 2):
+        temp_dir = tmp_path_factory.mktemp(f"temp_lhc_rdt_b{beam}", numbered=False)
+        linfile_dir = temp_dir / "lin_files"
+        model_dir = get_model_dir(beam=beam)
+        tmp_model_dir = temp_dir / model_dir.name
+
+        shutil.copytree(model_dir, tmp_model_dir)
+        decompress_model(tmp_model_dir)
+        run_harpy(beam=beam, linfile_dir=linfile_dir)
+
+        dfs[beam] = get_rdts_from_optics_analysis(
+            beam, linfile_dir=linfile_dir, output_dir=temp_dir, model_dir=tmp_model_dir
+        )
+    return dfs
+
+
+AMPLITUDE_TOLERANCES = {
+    1: {  # Beam 1: 0.6% and 4% for sextupoles and octupoles respectively
+        "sextupole": 6e-3,
+        "octupole": 4e-2,
+    },
+    2: {  # Beam 2: 0.6% and 5% for sextupoles and octupoles respectively
+        "sextupole": 6e-3,
+        "octupole": 5e-2,
+    },
+}
+PHASE_TOLERANCES = {
+    1: {  # Beam 1: 0.04 rad and 0.03 rad for orders 2 and 3 respectively
+        "sextupole": 4e-2,
+        "octupole": 3e-2,
+    },
+    2: {  # Beam 2: 0.07 rad and 0.05 rad for orders 2 and 3 respectively
+        "sextupole": 7e-2,
+        "octupole": 5e-2,
+    },
+}
+
+
+@pytest.mark.parametrize("rdt", get_rdt_names())
+@pytest.mark.parametrize("beam", [1, 2], ids=lambda x: f"Beam{x}")
+def test_lhc_rdts(beam: int, rdt: str, rdts_from_optics_analysis):
+    """Test the RDTs calculated by OMC3 against the analytical model and MAD-NG.
+
+    The test forces harpy to run for the selected RDTs to ensure the use of the
+    opposite_direction flag is consistent with the optics=True, beam=2 situation.
+
+    The configuration for the test is set in the file tests/inputs/lhc_rdts/rdt_constants.py:
+    - NTURNS = 1000
+    - KICK_AMP = 1e-3
+    - SEXTUPOLE_STRENGTH = 3e-5
+    - OCTUPOLE_STRENGTH = 3e-3
+
+    If this configuration is changed, the test needs to be updated accordingly by running
+    the create_data.py script in the same directory.
+
+    This configuration was specifically chosen to ensure that OMC3 is given a fair chance
+    to produce correct results. The main issues that can arise is detuning. The octupole
+    and sextupole strengths are chosen to be small enough to avoid large detuning effects,
+    but also large enough to produce a measurable effect. Furthermore, the larger each
+    strength, the more likely the regime will become nonlinear, which can lead to larger
+    errors in the RDT calculations.
+    """
+    # Retrieve the current RDT configuration
+    _, order = get_rdt_type(rdt)
+
+    # Now retrieve the optics calculations from the OMC3 analysis
+    rdt_dfs = rdts_from_optics_analysis[beam]
+    omc_df = rdt_dfs[rdt]  # get the result of a specific RDT
+
+    # Retrieve the MAD-NG results for the set of RDTs
+    file_suffix = get_file_suffix(beam)
+    ng_df = tfs.read(DATA_DIR / f"{MODEL_NG_PREFIX}_{file_suffix}.tfs", index="NAME")
+    ng_rdt = rdt.split("_")[0].upper()
+
+    # Reconstruct the complex numbers from the real and imaginary parts
+    ng_complex = ng_df[ng_rdt]
+    omc_complex = omc_df["REAL"] + 1j * omc_df["IMAG"]
+
+    # Calculate the OMC3 and MAD-NG amplitudes from the complex numbers
+    ng_amplitude = np.abs(ng_complex)
+    omc_amplitude = np.abs(omc_complex)
+
+    # Compare omc3 vs MAD-NG amplitudes, relative when above 1, absolute when below
+    ng_diff = ng_amplitude - omc_amplitude
+    ng_diff[ng_amplitude.abs() > 1] = (
+        ng_diff[ng_amplitude.abs() > 1] / ng_amplitude[ng_amplitude.abs() > 1]
+    )
+    assert ng_diff.abs().max() < AMPLITUDE_TOLERANCES[beam][order]
+
+    # Compare omc3 vs MAD-NG phases
+    ng_phase_diff = np.angle(ng_complex / omc_complex)
+    assert ng_phase_diff.max() < PHASE_TOLERANCES[beam][order]
+
+    # Check amplitude and phase vs real and imaginary calculations in omc3 is correct
+    omc_phase = np.angle(omc_complex) / PI2
+    diff = omc_phase - omc_df["PHASE"]
+    diff = (diff + 0.5) % 1 - 0.5
+    assert diff.abs().max() < 1e-11
+    assert np.allclose(omc_amplitude, omc_df["AMP"], rtol=1e-10)
+
+    # Now check the analytical model for sextupoles
+    # Analytical seems to disagree with MAD-NG and OMC3 for octupoles
+    if order == "sextupole":
+        analytical_df = tfs.read(
+            DATA_DIR / f"{MODEL_ANALYTICAL_PREFIX}_{file_suffix}.tfs", index="NAME"
+        )
+        analytical_complex = analytical_df[ng_rdt]
+
+        # Calculate the amplitudes
+        analytical_amplitude = np.abs(analytical_complex)
+        analytical_diff = analytical_amplitude - omc_amplitude
+
+        # Compare the amplitudes, relative when above 1, absolute when below
+        gt1 = analytical_amplitude.abs() > 1
+        analytical_diff[gt1] = analytical_diff[gt1] / analytical_amplitude[gt1]
+        assert analytical_diff.abs().max() < 1.1e-2
+
+        # Compare the phases
+        analytical_phase_diff = np.angle(analytical_complex / omc_complex)
+        assert analytical_phase_diff.max() < 6e-2