diff --git a/README.md b/README.md
index 283a333c47cd304089ffab416d8c6713a5b1f960..5eb6c63c46d60173605103b8108b9b11dd44c4c4 100644
--- a/README.md
+++ b/README.md
@@ -1,7 +1,9 @@
 mbtrack2
 ========
 
-[![Documentation Status](https://readthedocs.org/projects/mbtrack2/badge/?version=latest)](https://mbtrack2.readthedocs.io/en/latest/?badge=latest)
+![GitLab Release](https://img.shields.io/gitlab/v/release/184?gitlab_url=https%3A%2F%2Fgitlab.synchrotron-soleil.fr%2F)
+![PyPI - License](https://img.shields.io/pypi/l/mbtrack2)
+![Read the Docs](https://img.shields.io/readthedocs/mbtrack2)
 
 mbtrack2 is a coherent object-oriented framework written in python to work on collective effects in synchrotrons.
 
@@ -58,14 +60,22 @@ To test your installation run:
 from mbtrack2 import *
 ```
 
+### Using docker
+
+A docker image is available:
+
+```
+docker pull gitlab-registry.synchrotron-soleil.fr/pa/collective-effects/mbtrack2
+```
+
 References
 ----------
 If used in a publication, please cite mbtrack2 paper and the zenodo archive for the corresponding code version (and any other paper in this list for more specific features).
 
-[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.10871040.svg)](https://doi.org/10.5281/zenodo.10871040)
+[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.12749990.svg)](https://doi.org/10.5281/zenodo.12749990)
 
 ### mbtrack2 general features
 A. Gamelin, W. Foosang, and R. Nagaoka, “mbtrack2, a Collective Effect Library in Python”, presented at the 12th Int. Particle Accelerator Conf. (IPAC'21), Campinas, Brazil, May 2021, paper MOPAB070.
 
-### RF cavities with beam loading and RF feedbacks 
-Yamamoto, Naoto, Alexis Gamelin, and Ryutaro Nagaoka. "Investigation of Longitudinal Beam Dynamics With Harmonic Cavities by Using the Code Mbtrack." Proc. 10th International Partile Accelerator Conference (IPAC’19), Melbourne, Australia, 19-24 May 2019. 2019.
\ No newline at end of file
+### RF cavities with beam loading and RF feedbacks
+Yamamoto, Naoto, Alexis Gamelin, and Ryutaro Nagaoka. "Investigation of Longitudinal Beam Dynamics With Harmonic Cavities by Using the Code Mbtrack." Proc. 10th International Partile Accelerator Conference (IPAC’19), Melbourne, Australia, 19-24 May 2019. 2019.
diff --git a/mbtrack2/__init__.py b/mbtrack2/__init__.py
index 5ae28316144f7a09648c5da08b692ec81c3c3270..e478887b2fea2eeabdfbb44b0e5327efc3ade1c6 100644
--- a/mbtrack2/__init__.py
+++ b/mbtrack2/__init__.py
@@ -1,5 +1,5 @@
 # -*- coding: utf-8 -*-
-__version__ = "0.6.0"
+__version__ = "0.7.0"
 from mbtrack2.impedance import *
 from mbtrack2.instability import *
 from mbtrack2.tracking import *
diff --git a/mbtrack2/impedance/wakefield.py b/mbtrack2/impedance/wakefield.py
index 3325ed5a7060b5be0f69c6fc5024507c849e4da1..9b751fd20129155f8d370a6feb17f0e91071dd7e 100644
--- a/mbtrack2/impedance/wakefield.py
+++ b/mbtrack2/impedance/wakefield.py
@@ -1,6 +1,6 @@
 # -*- coding: utf-8 -*-
 """
-This module defines general classes to describes wakefields, impedances and 
+This module defines general classes to describes wakefields, impedances and
 wake functions.
 """
 
@@ -18,7 +18,7 @@ from scipy.interpolate import interp1d
 
 class ComplexData:
     """
-    Define a general data structure for a complex function based on a pandas 
+    Define a general data structure for a complex function based on a pandas
     DataFrame.
 
     Parameters
@@ -49,29 +49,29 @@ class ComplexData:
         """
         Method to add two structures. If the data don't have the same length,
         different cases are handled.
-        
+
         Parameters
         ----------
         structure_to_add : ComplexData object, int, float or complex.
             structure to be added.
         method : str, optional
-            manage how the addition is done, possibilties are: 
+            manage how the addition is done, possibilties are:
             -common: the common range of the index (variable) from the two
                 structures is used. The input data are cross-interpolated
                 so that the result data contains the points at all initial
                 points in the common domain.
             -extrapolate: extrapolate the value of both ComplexData.
-            -zero: outside of the common range, the missing values are zero. In 
+            -zero: outside of the common range, the missing values are zero. In
                 the common range, behave as "common" method.
         interp_kind : str, optional
             interpolation method which is passed to pandas and to
             scipy.interpolate.interp1d.
         index_name : str, optional
             name of the dataframe index passed to the method
-            
+
         Returns
         -------
-        ComplexData 
+        ComplexData
             Contains the sum of the two inputs.
         """
 
@@ -178,17 +178,17 @@ class ComplexData:
         """
         Interpolation of the data by calling the class to have a function-like
         behaviour.
-        
+
         Parameters
         ----------
         values : list or numpy array of complex, int or float
             values to be interpolated.
         interp_kind : str, optional
             interpolation method which is passed to scipy.interpolate.interp1d.
-            
+
         Returns
         -------
-        numpy array 
+        numpy array
             Contains the interpolated data.
         """
         real_func = interp1d(x=self.data.index,
@@ -338,7 +338,7 @@ class ComplexData:
 class WakeFunction(ComplexData):
     """
     Define a WakeFunction object based on a ComplexData object.
-    
+
     Parameters
     ----------
     variable : array-like
@@ -347,12 +347,12 @@ class WakeFunction(ComplexData):
         Wake function values in [V/C] for "long" and in [V/C/m] for others.
     component_type : str, optinal
         Type of the wake function: "long", "xdip", "xquad", ...
-    
+
     Attributes
     ----------
     data : DataFrame
     wake_type : str
-    
+
     Methods
     -------
     to_impedance()
@@ -432,9 +432,9 @@ class WakeFunction(ComplexData):
         """
         Return an Impedance object from the WakeFunction data.
         The WakeFunction data is assumed to be sampled equally.
-        
+
         If both sigma and mu are given, the impedance data is divided by the
-        Fourier transform of a Gaussian distribution. It then can be used to go 
+        Fourier transform of a Gaussian distribution. It then can be used to go
         from wake potential data to impedance.
 
         Parameters
@@ -496,11 +496,11 @@ class WakeFunction(ComplexData):
     def deconvolution(self, freq_lim, sigma, mu, nout=None):
         """
         Compute a wake function from wake potential data.
-                
-        The present data loaded into the WakeFunction object is assumed to be 
-        an uniformly sampled wake potential from a Gaussian bunch distribution 
+
+        The present data loaded into the WakeFunction object is assumed to be
+        an uniformly sampled wake potential from a Gaussian bunch distribution
         with RMS bunch length sigma and offset mu.
-        
+
         The offset mu depends on the code used to compute the wake potential:
             - for CST, it is the first time sample, so self.data.index[0].
             - for ABCI, it is "-1 * ISIG * SIG / c".
@@ -546,25 +546,25 @@ class WakeFunction(ComplexData):
 
     def loss_factor(self, sigma):
         """
-        Compute the loss factor or the kick factor for a Gaussian bunch. 
+        Compute the loss factor or the kick factor for a Gaussian bunch.
         Formulas from Eq. (5.6), Eq. (5.12) and Eq. (5.17) of [1].
-        
+
         Parameters
         ----------
         sigma : float
             RMS bunch length in [s]
-        
+
         Returns
         -------
         kloss: float
-            Loss factor in [V/C] or kick factor in [V/C/m] depanding on the 
+            Loss factor in [V/C] or kick factor in [V/C/m] depanding on the
             impedance type.
-            
+
         References
         ----------
-        [1] : Zotter, Bruno W., and Semyon A. Kheifets. Impedances and wakes 
+        [1] : Zotter, Bruno W., and Semyon A. Kheifets. Impedances and wakes
         in high-energy particle accelerators. World Scientific, 1998.
-        
+
         """
         time = self.data.index
         S = 1 / (2 * np.sqrt(np.pi) * sigma) * np.exp(-1 * time**2 /
@@ -577,7 +577,7 @@ class WakeFunction(ComplexData):
 class Impedance(ComplexData):
     """
     Define an Impedance object based on a ComplexData object.
-    
+
     Parameters
     ----------
     variable : array-like
@@ -586,12 +586,12 @@ class Impedance(ComplexData):
         Impedance values in [Ohm].
     component_type : str, optinal
         Type of the impedance: "long", "xdip", "xquad", ...
-    
+
     Attributes
     ----------
     data : DataFrame
     impedance_type : str
-    
+
     Methods
     -------
     loss_factor(sigma)
@@ -600,8 +600,6 @@ class Impedance(ComplexData):
         Return a WakeFunction object from the impedance data.
     plot()
         Plot the impedance data.
-    to_pycolleff()
-        Convenience method to export impedance to pycolleff.
     """
 
     def __init__(self,
@@ -673,20 +671,20 @@ class Impedance(ComplexData):
 
     def loss_factor(self, sigma):
         """
-        Compute the loss factor or the kick factor for a Gaussian bunch. 
+        Compute the loss factor or the kick factor for a Gaussian bunch.
         Formulas from Eq. (2) p239 and Eq.(7) p240 of [1].
-        
+
         Parameters
         ----------
         sigma : float
             RMS bunch length in [s]
-        
+
         Returns
         -------
         kloss: float
-            Loss factor in [V/C] or kick factor in [V/C/m] depanding on the 
+            Loss factor in [V/C] or kick factor in [V/C/m] depanding on the
             impedance type.
-            
+
         References
         ----------
         [1] : Handbook of accelerator physics and engineering, 3rd printing.
@@ -718,7 +716,7 @@ class Impedance(ComplexData):
         """
         Return a WakeFunction object from the impedance data.
         The impedance data is assumed to be sampled equally.
-    
+
         Parameters
         ----------
         nout : int or float, optional
@@ -727,11 +725,11 @@ class Impedance(ComplexData):
             input is padded with zeros. If nout < n, the input it cropped.
             Note that, for any nout value, nout//2+1 input points are required.
         trim : bool or float, optional
-            If True, the pseudo wake function is trimmed at time=0 and is forced 
-            to zero where time<0. 
+            If True, the pseudo wake function is trimmed at time=0 and is forced
+            to zero where time<0.
             If False, the original result is preserved.
-            If a float is given, the pseudo wake function is trimmed from 
-            time <= trim to 0. 
+            If a float is given, the pseudo wake function is trimmed from
+            time <= trim to 0.
         """
 
         Z0 = (self.data['real'] + self.data['imag'] * 1j)
@@ -783,49 +781,19 @@ class Impedance(ComplexData):
         ax.set_ylabel(label)
         return ax
 
-    def to_pycolleff(self):
-        """
-        Convenience method to export impedance to pycolleff.
-        Only implemented for longitudinal impedance.
-
-        Returns
-        -------
-        imp : pycolleff.longitudinal_equilibrium.ImpedanceSource
-            pycolleff ImpedanceSource object
-        """
-        # pycolleff is a soft dependency
-        from pycolleff.longitudinal_equilibrium import ImpedanceSource
-
-        # avoid circular import
-        from mbtrack2.utilities.misc import double_sided_impedance
-        imp = ImpedanceSource()
-        if self.component_type == "long":
-            double_sided_impedance(self)
-            # Negative imag sign convention !
-            imp.zl_table = self.data["real"].to_numpy(
-            ) - 1j * self.data["imag"].to_numpy()
-            imp.ang_freq_table = self.data.index.to_numpy() * 2 * np.pi
-        else:
-            raise NotImplementedError()
-
-        imp.calc_method = ImpedanceSource.Methods.ImpedanceDFT
-        imp.active_passive = ImpedanceSource.ActivePassive.Passive
-
-        return imp
-
 
 class WakeField:
     """
-    Defines a WakeField which corresponds to a single physical element which 
-    produces different types of wakes, represented by their Impedance or 
+    Defines a WakeField which corresponds to a single physical element which
+    produces different types of wakes, represented by their Impedance or
     WakeFunction objects.
-    
+
     Parameters
     ----------
     structure_list : list, optional
         list of Impedance/WakeFunction objects to add to the WakeField.
     name : str, optional
-    
+
     Attributes
     ----------
     impedance_components : array of str
@@ -834,7 +802,7 @@ class WakeField:
         WakeFunction components present for this element.
     components : array of str
         Impedance or WakeFunction components present for this element.
-        
+
     Methods
     -------
     append_to_model(structure_to_add)
@@ -984,7 +952,7 @@ class WakeField:
     def save(self, file):
         """
         Save WakeField element to file.
-        
+
         The same pandas version is needed on both saving and loading computer
         for the pickle to work.
 
@@ -1025,7 +993,7 @@ class WakeField:
     def load(file):
         """
         Load WakeField element from file.
-        
+
         The same pandas version is needed on both saving and loading computer
         for the pickle to work.
 
@@ -1088,7 +1056,7 @@ class WakeField:
     def add_several_wakefields(wakefields, beta):
         """
         Add a list of WakeFields taking into account beta functions.
-        
+
         Parameters
         ----------
         wakefields : list of WakeField
diff --git a/mbtrack2/tracking/rf.py b/mbtrack2/tracking/rf.py
index e39271204bc65fd8f015e43a34f4a21064e032d6..4de4aa28530c1f3942fa81d47a11f2fb19670ffa 100644
--- a/mbtrack2/tracking/rf.py
+++ b/mbtrack2/tracking/rf.py
@@ -1,6 +1,6 @@
 # -*- coding: utf-8 -*-
 """
-This module handles radio-frequency (RF) cavitiy elements. 
+This module handles radio-frequency (RF) cavitiy elements.
 """
 
 from functools import partial
@@ -18,7 +18,7 @@ class RFCavity(Element):
     """
     Perfect RF cavity class for main and harmonic RF cavities.
     Use cosine definition.
-    
+
     Parameters
     ----------
     ring : Synchrotron object
@@ -42,7 +42,7 @@ class RFCavity(Element):
         Tracking method for the element.
         No bunch to bunch interaction, so written for Bunch objects and
         @Element.parallel is used to handle Beam objects.
-        
+
         Parameters
         ----------
         bunch : Bunch or Beam object
@@ -58,15 +58,15 @@ class RFCavity(Element):
 class CavityResonator():
     """Cavity resonator class for active or passive RF cavity with beam
     loading or HOM, based on [1,2].
-    
+
     Use cosine definition.
-    
-    If used with mpi, beam.mpi.share_distributions must be called before the 
+
+    If used with mpi, beam.mpi.share_distributions must be called before the
     track method call.
-    
+
     Different kind of RF feeback and loops can be added using:
         cavity_resonator.feedback.append(loop)
-    
+
     Parameters
     ----------
     ring : Synchrotron object
@@ -89,11 +89,11 @@ class CavityResonator():
     theta : float, optional
         Total cavity phase objective value in [rad].
     n_bin : int, optional
-        Number of bins used for the beam loading computation. 
-        Only used if MPI is not used, otherwise n_bin must be specified in the 
+        Number of bins used for the beam loading computation.
+        Only used if MPI is not used, otherwise n_bin must be specified in the
         beam.mpi.share_distributions method.
         The default is 75.
-        
+
     Attributes
     ----------
     beam_phasor : complex
@@ -148,7 +148,7 @@ class CavityResonator():
     distance : array
         Distance between bunches.
     valid_bunch_index : array
-    
+
     Methods
     -------
     Vbr(I0)
@@ -193,17 +193,15 @@ class CavityResonator():
         Force feedback update from current CavityResonator parameters.
     sample_voltage()
         Sample the voltage seen by a zero charge particle during an RF period.
-    to_pycolleff()
-        Convenience method to export a CavityResonator to pycolleff.
-    
+
     References
     ----------
-    [1] Wilson, P. B. (1994). Fundamental-mode rf design in e+ e− storage ring 
-    factories. In Frontiers of Particle Beams: Factories with e+ e-Rings 
+    [1] Wilson, P. B. (1994). Fundamental-mode rf design in e+ e− storage ring
+    factories. In Frontiers of Particle Beams: Factories with e+ e-Rings
     (pp. 293-311). Springer, Berlin, Heidelberg.
-    
-    [2] Naoto Yamamoto, Alexis Gamelin, and Ryutaro Nagaoka. "Investigation 
-    of Longitudinal Beam Dynamics With Harmonic Cavities by Using the Code 
+
+    [2] Naoto Yamamoto, Alexis Gamelin, and Ryutaro Nagaoka. "Investigation
+    of Longitudinal Beam Dynamics With Harmonic Cavities by Using the Code
     Mbtrack." IPAC’19, Melbourne, Australia, 2019.
     """
 
@@ -262,11 +260,11 @@ class CavityResonator():
     def track(self, beam):
         """
         Track a Beam object through the CavityResonator object.
-        
+
         Can be used with or without mpi.
         If used with mpi, beam.mpi.share_distributions must be called before.
-        
-        The beam phasor is given at t=0 (synchronous particle) of the first 
+
+        The beam phasor is given at t=0 (synchronous particle) of the first
         non empty bunch.
 
         Parameters
@@ -366,8 +364,8 @@ class CavityResonator():
         """
         Initialize the beam phasor for a given beam distribution using a
         tracking like method.
-        
-        Follow the same steps as the track method but in the "rf" reference 
+
+        Follow the same steps as the track method but in the "rf" reference
         frame and without any modifications on the beam.
 
         Parameters
@@ -426,7 +424,7 @@ class CavityResonator():
 
     def phasor_decay(self, time, ref_frame="beam"):
         """
-        Compute the beam phasor decay during a given time span, assuming that 
+        Compute the beam phasor decay during a given time span, assuming that
         no particles are crossing the cavity during the time span.
 
         Parameters
@@ -450,9 +448,9 @@ class CavityResonator():
                     charge_per_mp,
                     ref_frame="beam"):
         """
-        Compute the beam phasor evolution during the crossing of a bunch using 
+        Compute the beam phasor evolution during the crossing of a bunch using
         an analytic formula [1].
-        
+
         Assume that the phasor decay happens before the beam loading.
 
         Parameters
@@ -465,11 +463,11 @@ class CavityResonator():
             Charge per macro-particle in [C].
         ref_frame : string, optional
             Reference frame to be used, can be "beam" or "rf".
-            
+
         References
         ----------
         [1] mbtrack2 manual.
-            
+
         """
         if ref_frame == "beam":
             delta = self.wr
@@ -494,13 +492,13 @@ class CavityResonator():
         """
         Initialize the beam phasor for a given beam distribution using an
         analytic formula [1].
-        
+
         No modifications on the Beam object.
 
         Parameters
         ----------
         beam : Beam object
-            
+
         References
         ----------
         [1] mbtrack2 manual.
@@ -640,7 +638,7 @@ class CavityResonator():
 
     @property
     def Rs_per_cavity(self):
-        """Shunt impedance of a single cavity in [Ohm], defined as 
+        """Shunt impedance of a single cavity in [Ohm], defined as
         0.5*Vc*Vc/Pc."""
         return self._Rs_per_cavity
 
@@ -844,7 +842,7 @@ class CavityResonator():
 
     def set_generator(self, I0):
         """
-        Set generator parameters (Pg, Vgr, theta_gr, Vg and theta_g) for a 
+        Set generator parameters (Pg, Vgr, theta_gr, Vg and theta_g) for a
         given current and set of parameters.
 
         Parameters
@@ -881,14 +879,14 @@ class CavityResonator():
 
     def plot_phasor(self, I0):
         """
-        Plot phasor diagram showing the vector addition of generator and beam 
+        Plot phasor diagram showing the vector addition of generator and beam
         loading voltage.
 
         Parameters
         ----------
         I0 : float
             Beam current in [A].
-            
+
         Returns
         -------
         Figure.
@@ -952,12 +950,12 @@ class CavityResonator():
                       modes=None,
                       bool_return=False):
         """
-        Check Coupled-Bunch-Instability stability, 
+        Check Coupled-Bunch-Instability stability,
         Effect of Direct RF feedback is not included.
 
-        This method is a wraper around lcbi_growth_rate to caluclate the CBI 
+        This method is a wraper around lcbi_growth_rate to caluclate the CBI
         growth rate from the CavityResonator own impedance.
-        
+
         See lcbi_growth_rate for details.
 
         Parameters
@@ -966,7 +964,7 @@ class CavityResonator():
             Beam current in [A].
         synchrotron_tune : float, optinal
             Fractional number of longitudinal tune.
-            If None, synchrotron_tune is computed using self.Vc as total 
+            If None, synchrotron_tune is computed using self.Vc as total
             voltage.
             Default is None.
         modes : float or list of float, optional
@@ -975,17 +973,17 @@ class CavityResonator():
             Default is None.
         bool_return : bool, optional
             If True:
-                - return True if the most unstable mode is stronger than 
+                - return True if the most unstable mode is stronger than
                 longitudinal damping time.
                 - return False if it is not the case.
             Default is False.
-            
+
         Returns
         -------
         growth_rate : float
             Maximum coupled bunch instability growth rate in [s-1].
         mu : int
-            Coupled bunch mode number corresponding to the maximum coupled bunch 
+            Coupled bunch mode number corresponding to the maximum coupled bunch
             instability growth rate.
 
         """
@@ -1105,7 +1103,7 @@ class CavityResonator():
         n_points : int or float, optional
             Number of sample points. The default is 1e4.
         index : int, optional
-            RF bucket number to sample. Be carful if index > 0 as no new beam 
+            RF bucket number to sample. Be carful if index > 0 as no new beam
             loading is not taken into account here.
             The default is 0.
 
@@ -1146,54 +1144,13 @@ class CavityResonator():
 
         return pos, voltage_rec
 
-    def to_pycolleff(self, Impedance=True):
-        """
-        Convenience method to export a CavityResonator to pycolleff.
-        
-        Parameters
-        ----------
-        Impedance : bool, optional
-            If True, export as impedance (i.e. ImpedanceSource.Methods.ImpedanceDFT).
-            If False, export as wake (i.e. ImpedanceSource.Methods.Wake).
-            Default is True.
-
-        Returns
-        -------
-        cav : pycolleff ImpedanceSource object
-
-        """
-        # pycolleff is a soft dependency
-        from pycolleff.longitudinal_equilibrium import ImpedanceSource
-
-        cav = ImpedanceSource()
-        cav.harm_rf = self.m
-        cav.Q = self.QL
-        RoverQ = self.RL / self.QL
-        cav.shunt_impedance = RoverQ * cav.Q
-        cav.ang_freq_rf = self.ring.omega1
-        cav.ang_freq = cav.harm_rf * cav.ang_freq_rf
-        cav.detune_w = 2 * np.pi * self.detune
-        if self.Vg != 0:
-            cav.active_passive = ImpedanceSource.ActivePassive.Active
-            if Impedance:
-                raise NotImplementedError()
-            else:
-                cav.calc_method = ImpedanceSource.Methods.Wake
-        else:
-            cav.active_passive = ImpedanceSource.ActivePassive.Passive
-            if Impedance:
-                cav.calc_method = ImpedanceSource.Methods.ImpedanceDFT
-            else:
-                cav.calc_method = ImpedanceSource.Methods.Wake
-        return cav
-
 
 class ProportionalLoop():
     """
     Proportional feedback loop to control a CavityResonator amplitude and phase.
-    
+
     Feedback setpoints are cav_res.Vc and cav_res.theta.
-    
+
     The loop must be added to the CavityResonator object using:
         cav_res.feedback.append(loop)
 
@@ -1248,9 +1205,9 @@ class TunerLoop():
     """
     Cavity tuner loop used to control a CavityResonator tuning (psi or detune)
     in order to keep the phase between cavity and generator current constant.
-    
+
     Only a proportional controller is assumed.
-    
+
     The loop must be added to the CavityResonator object using:
         cav_res.feedback.append(loop)
 
@@ -1263,7 +1220,7 @@ class TunerLoop():
         Proportional gain of the tuner loop.
         If not specified, 0.01 is used.
     avering_period:
-        Period during which the phase difference is monitored and averaged. 
+        Period during which the phase difference is monitored and averaged.
         Then the feedback correction is applied every avering_period turn.
         Unit is turn number.
         A value longer than one synchrotron period (1/fs) is recommended.
@@ -1271,7 +1228,7 @@ class TunerLoop():
         too fast compared to the actual situation.
     offset : float
         Tuning offset in [rad].
-        
+
     """
 
     def __init__(self, ring, cav_res, gain=0.01, avering_period=0, offset=0):
@@ -1291,11 +1248,11 @@ class TunerLoop():
     def track(self):
         """
         Tracking method for the tuner loop.
-    
+
         Returns
         -------
         None.
-    
+
         """
         if self.count == self.avering_period:
             diff = self.diff / self.avering_period - self.offset
@@ -1309,35 +1266,35 @@ class TunerLoop():
 
 class ProportionalIntegralLoop():
     """
-    Proportional Integral (PI) loop to control a CavityResonator amplitude and 
+    Proportional Integral (PI) loop to control a CavityResonator amplitude and
     phase via generator current (ig) [1,2].
-    
+
     Feedback reference targets (setpoints) are cav_res.Vc and cav_res.theta.
-    
-    The ProportionalIntegralLoop should be initialized only after generator 
+
+    The ProportionalIntegralLoop should be initialized only after generator
     parameters are set.
-    
+
     The loop must be added to the CavityResonator object using:
         cav_res.feedback.append(loop)
-        
-    When the reference target is changed, it might be needed to re-initialize 
+
+    When the reference target is changed, it might be needed to re-initialize
     the feedforward constant by calling init_FFconst().
 
     Parameters
     ----------
     ring : Synchrotron object
     cav_res : CavityResonator object
-        CavityResonator in which the loop will be added. 
+        CavityResonator in which the loop will be added.
     gain : float list
-        Proportional gain (Pgain) and integral gain (Igain) of the feedback 
+        Proportional gain (Pgain) and integral gain (Igain) of the feedback
         system in the form of a list [Pgain, Igain].
-        In case of normal conducting cavity (QL~1e4), the Pgain of ~1.0 and 
+        In case of normal conducting cavity (QL~1e4), the Pgain of ~1.0 and
         Igain of ~1e4(5) are usually used.
-        In case of super conducting cavity (QL > 1e6), the Pgain of ~100 
+        In case of super conducting cavity (QL > 1e6), the Pgain of ~100
         can be used.
-        In a "bad" parameter set, unstable oscillation of the cavity voltage 
+        In a "bad" parameter set, unstable oscillation of the cavity voltage
         can be caused. So, a parameter scan of the gain should be made.
-        Igain * ring.T1 / dtau is Ki defined as a the coefficients for 
+        Igain * ring.T1 / dtau is Ki defined as a the coefficients for
         integral part in of [1], where dtau is a clock period of the PI controller.
     sample_num : int
         Number of bunch over which the mean cavity voltage is computed.
@@ -1347,7 +1304,7 @@ class ProportionalIntegralLoop():
         Time interval between two cavity voltage monitoring and feedback.
         Units are in bucket numbers.
     delay : int
-        Loop delay of the PI feedback system. 
+        Loop delay of the PI feedback system.
         Units are in bucket numbers.
     IIR_cutoff : float, optinal
         Cutoff frequency of the IIR filter in [Hz].
@@ -1355,12 +1312,12 @@ class ProportionalIntegralLoop():
         Default is 0.
     FF : bool, optional
         Boolean switch to use feedforward constant.
-        True is recommended to prevent a cavity voltage drop in the beginning 
+        True is recommended to prevent a cavity voltage drop in the beginning
         of the tracking.
-        In case of small Pgain (QL ~ 1e4), cavity voltage drop may cause 
+        In case of small Pgain (QL ~ 1e4), cavity voltage drop may cause
         beam loss due to static Robinson.
         Default is True.
-        
+
     Methods
     -------
     track()
@@ -1385,7 +1342,7 @@ class ProportionalIntegralLoop():
     Notes
     -----
     Assumption : delay >> every ~ sample_num
-    
+
     Adjusting ig(Vg) parameter to keep the cavity voltage constant (to target values).
     The "track" method is
        1) Monitoring the cavity voltage phasor
@@ -1409,14 +1366,14 @@ class ProportionalIntegralLoop():
                      is one reasonable parameter set.
         The practical clock period is 13ns.
             ==> Igain_PF = Igain_mbtrack * Trf / 13ns = Igain_mbtrack * 0.153
-     
+
     References
     ----------
     [1] : https://en.wikipedia.org/wiki/PID_controller
-    [2] : Yamamoto, N., Takahashi, T., & Sakanaka, S. (2018). Reduction and 
-    compensation of the transient beam loading effect in a double rf system 
+    [2] : Yamamoto, N., Takahashi, T., & Sakanaka, S. (2018). Reduction and
+    compensation of the transient beam loading effect in a double rf system
     of synchrotron light sources. PRAB, 21(1), 012001.
-    
+
     """
 
     def __init__(self,
@@ -1506,8 +1463,8 @@ class ProportionalIntegralLoop():
     def init_Ig2Vg_matrix(self):
         """
         Initialize matrix for Ig2Vg_matrix.
-        
-        Shoud be called before first use of Ig2Vg_matrix and after each cavity 
+
+        Shoud be called before first use of Ig2Vg_matrix and after each cavity
         parameter change.
         """
         k = np.arange(0, self.ring.h)
@@ -1529,7 +1486,7 @@ class ProportionalIntegralLoop():
     def Ig2Vg_matrix(self):
         """
         Return Vg from Ig using matrix formalism.
-        Warning: self.init_Ig2Vg should be called after each CavityResonator 
+        Warning: self.init_Ig2Vg should be called after each CavityResonator
         parameter change.
         """
         generator_phasor_record = (
@@ -1541,7 +1498,7 @@ class ProportionalIntegralLoop():
     def Ig2Vg(self):
         """
         Go from Ig to Vg.
-        
+
         Apply new values to cav_res.generator_phasor_record, cav_res.Vg and
         cav_res.theta_g from ig_phasor_record.
         """
@@ -1567,7 +1524,7 @@ class ProportionalIntegralLoop():
         cutoff : float
             Cutoff frequency of the IIR filter in [Hz].
             If 0, cutoff frequency is infinity.
-            
+
         """
         if cutoff == 0:
             self.IIRcoef = 1.0
@@ -1599,18 +1556,18 @@ class ProportionalIntegralLoop():
 class DirectFeedback(ProportionalIntegralLoop):
     """
     Direct RF FB (DFB) via generator current (ig) using PI controller [1,2].
-    
-    The DirectFeedback inherits from ProportionalIntegralLoop so all 
+
+    The DirectFeedback inherits from ProportionalIntegralLoop so all
     mandatory parameters from ProportionalIntegralLoop should be passed as
     kwargs when initializing a DirectFeedback object.
-        
+
     To avoid cavity-beam unmatching (large synchrotron oscilation of beam),
-    the CavityResonator generator parameters should be set before 
+    the CavityResonator generator parameters should be set before
     initialization.
-    
+
     The loop must be added to the CavityResonator object using:
         cav_res.feedback.append(loop)
-    
+
     Parameters
     ----------
     DFB_gain : float
@@ -1633,7 +1590,7 @@ class DirectFeedback(ProportionalIntegralLoop):
         Units are in bucket numbers.
         If None, value from the ProportionalIntegralLoop is used.
         Default is None.
-        
+
     Attributes
     ----------
     phase_shift : float
@@ -1646,7 +1603,7 @@ class DirectFeedback(ProportionalIntegralLoop):
         Return the gain gamma of the DFB.
     DFB_Rs : float
         Return the shunt impedance of the DFB in [ohm].
-    
+
     Methods
     -------
     DFB_parameter_set(DFB_gain, DFB_phase_shift)
@@ -1657,14 +1614,14 @@ class DirectFeedback(ProportionalIntegralLoop):
         Return the generator voltage main and DFB components in [V].
     DFB_fs()
         Return the modified synchrotron frequency in [Hz].
-        
+
     References
     ----------
-    [1] : Akai, K. (2022). Stability analysis of rf accelerating mode with 
-    feedback loops under heavy beam loading in SuperKEKB. PRAB, 25(10), 
+    [1] : Akai, K. (2022). Stability analysis of rf accelerating mode with
+    feedback loops under heavy beam loading in SuperKEKB. PRAB, 25(10),
     102002.
-    [2] : N. Yamamoto et al. (2023). Stability survey of a double RF system 
-    with RF feedback loops for bunch lengthening in a low-emittance synchrotron 
+    [2] : N. Yamamoto et al. (2023). Stability survey of a double RF system
+    with RF feedback loops for bunch lengthening in a low-emittance synchrotron
     ring. In Proc. IPAC'23. doi:10.18429/JACoW-IPAC2023-WEPL161
 
     """
@@ -1731,7 +1688,7 @@ class DirectFeedback(ProportionalIntegralLoop):
     def DFB_psi(self):
         """
         Return detuning angle value with Direct RF feedback in [rad].
-    
+
         Fig.4 of ref [1].
         """
         return (np.angle(np.mean(self.cav_res.cavity_phasor_record)) -
@@ -1741,7 +1698,7 @@ class DirectFeedback(ProportionalIntegralLoop):
     def DFB_alpha(self):
         """
         Return the amplitude ratio alpha of the DFB.
-        
+
         Near Eq. (13) of [1].
         """
         fac = np.abs(
@@ -1752,7 +1709,7 @@ class DirectFeedback(ProportionalIntegralLoop):
     def DFB_gamma(self):
         """
         Return the gain gamma of the DFB.
-        
+
         Near Eq. (13) of [1].
         """
         fac = np.abs(
@@ -1763,7 +1720,7 @@ class DirectFeedback(ProportionalIntegralLoop):
     def DFB_Rs(self):
         """
         Return the shunt impedance of the DFB in [ohm].
-        
+
         Eq. (15) of [1].
         """
         return self.cav_res.Rs / (1 + self.DFB_gamma * np.cos(self.DFB_psi))
diff --git a/mbtrack2/tracking/synchrotron.py b/mbtrack2/tracking/synchrotron.py
index e9df8ffe23882de2af3d7631588201c77f9b70bc..15235ecd1b60f38a4984074ab71c87c1fba54a85 100644
--- a/mbtrack2/tracking/synchrotron.py
+++ b/mbtrack2/tracking/synchrotron.py
@@ -11,10 +11,10 @@ from scipy.constants import c, e
 class Synchrotron:
     """
     Synchrotron class to store main properties.
-    
-    Optional parameters are optional only if the Optics object passed to the 
+
+    Optional parameters are optional only if the Optics object passed to the
     class uses a loaded lattice.
-    
+
     Parameters
     ----------
     h : int
@@ -44,25 +44,25 @@ class Synchrotron:
     U0 : float, optional
         Energy loss per turn in [eV].
     adts : list of arrays or None, optional
-        List that contains arrays of polynomial's coefficients, in decreasing 
-        powers, used to determine Amplitude-Dependent Tune Shifts (ADTS). 
+        List that contains arrays of polynomial's coefficients, in decreasing
+        powers, used to determine Amplitude-Dependent Tune Shifts (ADTS).
         The order of the elements strictly needs to be
         [coef_xx, coef_yx, coef_xy, coef_yy], where x and y denote the horizontal
         and the vertical plane, respectively, and coef_PQ means the polynomial's
         coefficients of the ADTS in plane P due to the amplitude in plane Q.
 
-        For example, if the tune shift in y due to the Courant-Snyder invariant 
-        Jx is characterized by the equation dQy(x) = 3*Jx**2 + 2*Jx + 1, then 
+        For example, if the tune shift in y due to the Courant-Snyder invariant
+        Jx is characterized by the equation dQy(x) = 3*Jx**2 + 2*Jx + 1, then
         coef_yx takes the form np.array([3, 2, 1]).
-        
+
         Use None, to exclude the ADTS from calculation.
     mcf_order : array, optional
         Higher-orders momentum compaction factor in decreasing powers.
         Input here the coefficients considering only the derivatives and which
         do not include any factorial factor.
-        The 1st order should be included in the array and be at the last 
+        The 1st order should be included in the array and be at the last
         position.
-        
+
     Attributes
     ----------
     T0 : float
@@ -92,8 +92,8 @@ class Synchrotron:
     long_gamma : float
         Longitudinal gamma Twiss parameter at the tracking location in [s-1].
         Initialized at zero.
-        
-        
+
+
     Methods
     -------
     synchrotron_tune(V)
@@ -104,20 +104,19 @@ class Synchrotron:
         Momentum compaction taking into account higher orders if provided in
         mcf_order.
     mcf(delta)
-        Momentum compaction factor taking into account higher orders if 
+        Momentum compaction factor taking into account higher orders if
         provided in mcf_order.
     get_adts()
-        Compute and add Amplitude-Dependent Tune Shifts (ADTS) sextupolar 
+        Compute and add Amplitude-Dependent Tune Shifts (ADTS) sextupolar
         componenet from AT lattice.
     get_mcf_order()
         Compute momentum compaction factor up to 3rd order from AT lattice.
     get_longitudinal_twiss(V)
-        Compute the longitudinal Twiss parameters and the synchrotron tune for 
+        Compute the longitudinal Twiss parameters and the synchrotron tune for
         single or multi-harmonic RF systems.
     to_pyat(Vrf)
         Return a pyAT simple_ring element from the Synchrotron element data.
-    to_pycolleff(I0, Vrf, bunch_number)
-        Return a pycolleff Ring element from the Synchrotron element data.
+
     """
 
     def __init__(self, h, optics, particle, **kwargs):
@@ -307,13 +306,13 @@ class Synchrotron:
         Parameters
         ----------
         position : float or array, optional
-            Longitudinal position in [m] where the beam size is computed. 
+            Longitudinal position in [m] where the beam size is computed.
             If None, the local values are used.
 
         Returns
         -------
         sigma : array
-            RMS beam size in [m] at position location or at local positon if 
+            RMS beam size in [m] at position location or at local positon if
             position is None.
 
         """
@@ -354,7 +353,7 @@ class Synchrotron:
     def synchrotron_tune(self, V):
         """
         Compute the (unperturbed) synchrotron tune from main RF voltage.
-        
+
         Parameters
         ----------
         V : float
@@ -364,7 +363,7 @@ class Synchrotron:
         -------
         tuneS : float
             Synchrotron tune.
-            
+
         """
         Vsum = V * np.sin(np.arccos(self.U0 / V))
         phi = np.arccos(1 - self.eta(0) * np.pi * self.h / self.E0 * Vsum)
@@ -373,7 +372,7 @@ class Synchrotron:
 
     def get_adts(self):
         """
-        Compute and add Amplitude-Dependent Tune Shifts (ADTS) sextupolar 
+        Compute and add Amplitude-Dependent Tune Shifts (ADTS) sextupolar
         componenet from AT lattice.
         """
         import at
@@ -391,7 +390,7 @@ class Synchrotron:
     def get_mcf_order(self, add=True, show_fit=False):
         """
         Compute momentum compaction factor up to 3rd order from AT lattice.
-        
+
         Parameters
         ----------
         add : bool, optional
@@ -401,7 +400,7 @@ class Synchrotron:
         show_fit : bool, optional
             If True, show a plot with the data and fit.
             The default is False.
-            
+
         Return
         ------
         pvalue : array, optional
@@ -441,13 +440,13 @@ class Synchrotron:
 
     def get_longitudinal_twiss(self, V, phase=None, harmonics=None, add=True):
         """
-        Compute the longitudinal Twiss parameters and the synchrotron tune for 
+        Compute the longitudinal Twiss parameters and the synchrotron tune for
         single or multi-harmonic RF systems.
-        
+
         Hypthesis:
         - Use the linear approximation (tau ~ 0).
         - For higher haromics only, cos(phi) ~ 0 must be true.
-        
+
         Parameters
         ----------
         V : float or array-like of float
@@ -478,7 +477,7 @@ class Synchrotron:
         long_alpha : float
             Longitudinal alpha Twiss parameter at the tracking location.
         long_beta : float
-            Longitudinal beta Twiss parameter at the tracking location in [s].   
+            Longitudinal beta Twiss parameter at the tracking location in [s].
         long_gamma : float
             Longitudinal gamma Twiss parameter at the tracking location in [s-1].
         """
@@ -516,7 +515,7 @@ class Synchrotron:
     def to_pyat(self, Vrf, harmonic_number=None, TimeLag=False):
         """
         Return a pyAT simple_ring element from the Synchrotron element data.
-        
+
         See pyAT documentation for informations about simple_ring.
 
         Parameters
@@ -527,7 +526,7 @@ class Synchrotron:
         harmonic_number : float or array-like of float, optional
             Harmonic numbers of the RF cavities. The default is None.
         TimeLag : float or array-like of float, optional
-            Set the timelag of the cavities in pyAT definition. 
+            Set the timelag of the cavities in pyAT definition.
             The default is False.
 
         Returns
@@ -587,48 +586,3 @@ class Synchrotron:
                                      U0=self.U0,
                                      TimeLag=TimeLag)
         return at_simple_ring
-
-    def to_pycolleff(self, I0, Vrf, bunch_number):
-        """
-        Return a pycolleff Ring element from the Synchrotron element data.
-
-        Parameters
-        ----------
-        I0 : float
-            Beam current in [A].
-        Vrf : float
-            RF voltage in [V].
-        bunch_number : int
-            Total number of bunches filled.
-
-        Returns
-        -------
-        ring : pycolleff.colleff.Ring
-            pycolleff Ring object.
-
-        """
-        # pycolleff is a soft dependency
-        from pycolleff.colleff import Ring
-
-        ring = Ring()
-        ring.version = 'from_mbtrack2'
-        ring.rf_freq = self.f1
-        ring.mom_comp = self.ac  # momentum compaction factor
-        ring.energy = self.E0  # energy [eV]
-        ring.tuney = self.tune[1]  # vertical tune
-        ring.tunex = self.tune[0]  # horizontal tune
-        ring.chromx = self.chro[0]  # horizontal chromaticity
-        ring.chromy = self.chro[1]  # vertical chromaticity
-        ring.harm_num = self.h  # harmonic Number
-        ring.num_bun = bunch_number  # number of bunches filled
-        ring.total_current = I0  # total current [A]
-        ring.sync_tune = self.synchrotron_tune(Vrf)  # synchrotron tune
-        ring.espread = self.sigma_delta
-        ring.bunlen = self.sigma_0 * c  # [m]
-        ring.damptx = self.tau[0]  # [s]
-        ring.dampty = self.tau[1]  # [s]
-        ring.dampte = self.tau[2]  # [s]
-        ring.en_lost_rad = self.U0  # [eV]
-        ring.gap_voltage = Vrf  # [V]
-
-        return ring
diff --git a/pyproject.toml b/pyproject.toml
index 8ac63cbe5d1eda70cb26ce009661bcc3ed3d5f44..5e2b88465bd69d95670f721ee64a30dbef2df0ec 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "mbtrack2"
-version = "0.6.0"
+version = "0.7.0"
 description = "A coherent object-oriented framework to work on collective effects in synchrotrons."
 authors = ["Alexis Gamelin <alexis.gamelin@synchrotron-soleil.fr>"]
 license = "BSD-3-Clause"
@@ -31,8 +31,8 @@ build-backend = "poetry.core.masonry.api"
 [tool.yapf]
 based_on_style = "pep8"
 arithmetic_precedence_indication = true
-blank_line_before_nested_class_or_def = true 
+blank_line_before_nested_class_or_def = true
 
 [tool.isort]
 multi_line_output = 3
-include_trailing_comma = true
\ No newline at end of file
+include_trailing_comma = true