add lifetime and momentum aperture

pyatsoleiloutils/ltandmomap.py

+"""This file contains functions to calculate de momentum aperture
+and lifetime of a ring
+"""
+
+import numpy
+import scipy
+import at
+import h5py
+from .pyatsoleiloutils import savemomaphdf5
+from .pyatsoleiloutils import readmomaphdf5
+
+__all__ = ['gettouscheklt_SOLEIL_II','getmomap_SOLEIL_II']
+
+def getmomaphdf5_SOLEIL_II(ring, foutn='momap.hdf5'):
+    """
+    This function calculates the momentum aperture of a ring
+    and saves the result in a hdf5 file
+    """
+    dpb, s_indexes, s_array = getmomap_SOLEIL_II(ring)
+    savemomaphdf5(
+            foutn,
+            dpb,
+            s_indexes,
+            s_array)
+ 
+def getmomap_SOLEIL_II(ring):
+    """
+    This function calculates the momentum aperture of a ring
+    """
+    # check the s coordinate along the machine
+    s_list = at.get_s_pos(ring, range(len(ring)))
+    #print(f'{len(ring)=}')
+    # choose unique s
+    refpts = numpy.arange(len(ring))
+    mask = [ring[r].Length > 0.0 for r in refpts]
+    refpts = refpts[mask]
+
+    # get momentum aperture
+    s_indexes = refpts
+    s_array = at.get_s_pos(ring, s_indexes)
+    dp_boundary, _, _ = at.get_momentum_acceptance(ring,
+            resolution=0.5e-3,
+            amplitude=0.10,
+            nturns=512,
+            grid_mode=at.GridMode(1),
+            refpts=s_indexes,
+            use_mp=True,
+            verbose=True,
+            divider=2
+            )
+
+    return dp_boundary, s_indexes, s_array
+
+
+def gettouscheklt_SOLEIL_II(ring, beammode='TDR',finname=None):
+    """
+    This function calculates the Touschek lifetime of a ring
+    using a list of beam parameters
+    """
+    beam_par = at.envelope_parameters(ring)
+    ex = beam_par.emittances[0]
+    nbunches = 416
+    i_firstmode = 500e-3 #[A]
+
+    if beammode == 'TDR':
+       ey = 0.3*ex  # SOLEIL II TDR emittance
+    else:
+       # 50% of emmittance on each plane
+       ey = 0.5*ex
+       ex = ey
+
+    if finnanme == None:
+       momap_pn, s_indexes, s_array = getmomap_SOLEIL_II(ring)
+    else:
+       momap_pn, s_indexes, s_array = readmomaphdf5(momapfile)
+
+    # get touschek lifetime
+    lt_tou, _, _ = at.touschek.get_lifetime(ring,
+            ey,
+            i_firstmode/nbunches,
+            emitx=ex,
+            sigs=beam_par.sigma_l,
+            sigp=beam_par.sigma_e,
+            momap=momap_pn,
+            refpts=s_indexes
+            )
+
+    return lt_tou
+