Robust linear coupling correction with N-turn maps

062801-1 The linear one-turn map of a storage ring contains coupling information on which a correction algorithm can be based. In principal, the one-turn matrix can be fitted from turn-by-turn data of beam position monitors. However, the signal-to-noise ratio of the coupling information can be greatly enhanced by fitting maps for larger turn numbers N. Furthermore, by using a number of beam position monitors in a region with only small coupling sources, the determination of the N-turn map can be made robust against failures of individual beam position monitors, and the signal-to-noise ratio of the coupling information is further enhanced. With the so obtained N-turn maps an automated global coupling correction is possible without the need for a tune change. This is demonstrated for the Brookhaven Relativistic Heavy Ion Collider where the implementation of the algorithm allows a global coupling correction within a few seconds at injection.


I. INTRODUCTION
Linear coupling, when large enough, can make it difficult or impossible to set tunes to values close to the coupling resonance Q x Q y .These tunes are desirable since the resonance density in this area is low.When a tune feedback is used, linear coupling may prevent the feedback system from attaining the set values and cause the loops to open.A fast and robust method for decoupling is therefore desirable.
A widely used method to measure global linear coupling is to move the tunes until the minimum tune separation Q min jQ x ÿ Q y j min is reached [1].A coupling correction is then performed by scanning skew quadrupole corrector settings to minimize Q min .This approach has several disadvantages.First, due to its iterative nature the procedure is time consuming.This can result in time lost for luminosity production, or slow performance degradation if the correction is not performed routinely.Second, in a coupled machine scanning one tune may push the other tune close to a limiting resonance.In superconducting machines, the resulting beam loss can lead to an aborted store to protect the magnets from quenching.A method that does not change the tunes is therefore called for.Third, a tune and skew quadrupole scanning method cannot be practically applied to a machine during an energy ramp.The decoupling method presented here overcomes all these shortcomings.
In lepton machines, a coupling measurement and correction can be based on a continuous beam excitation with a shaker (see, for example, Refs.[2 -4]).However, the beam excitation with a shaker leads to unacceptable emittance growth in hadron machines.Apart from the tune scanning method mentioned above, other coupling measurements and correction algorithms were implemented in hadron machines [5,6].But in these cases, hardware beyond a turn-by-turn orbit acquisition system was used.
The measurement and correction of linear coupling in storage rings was analyzed by a number of authors.A short collection of this work is given as Refs.[1][2][3][4][5][6][7][8][9][10][11][12].The purpose of this article is not to repeat such an analysis.Rather, based on previous work, an algorithm is presented to measure and correct linear coupling globally from turn-by-turn orbit data.We assume here that the orbit data are taken after a transverse kick was applied.In principle the kick size can be small enough to avoid sizable emittance growth.Emphasis is put on robustness against failures of individual beam position monitors and random errors in monitor readings.Based on fitted N-turn maps the minimum tune approach Q min jQ x ÿ Q y j min , a common measure of linear coupling, and skew corrector settings can be obtained without a tune change.The algorithm lends itself to full automation and allows a coupling correction within seconds.When turn-by-turn data can be obtained during an energy ramp, skew corrector settings can be computed for the next ramp.
In the following we review the matrix description of linear coupling.The construction of turn-by-turn 4-vectors is shown, from which N-turn maps can be fitted.The minimum tune separation Q min is computed from the N-turn map, as well as skew corrector settings to minimize it.The implementation of an automatic coupling correction, based on this method, is presented for the Relativistic Heavy Ion Collider (RHIC).
The 4 4 matrix M can be written in terms of 2 2 matrices as The machine is said to be globally decoupled if B 0, which implies C 0. If the machine is globally decoupled at one observation point, it is not necessarily decoupled at another observation point.B 0 does not mean that the eigenmotions are in the horizontal and vertical planes, but that the resonant energy excitation between one eigenmode and the other is suppressed.Furthermore, if the linear coupling is caused by a number of small sources, rather than a few large ones, global decoupling at any observation point will usually lead to a machine that is almost globally decoupled at any other observation point [11].We denote by Q A ; Q D the eigentunes of the map (2), where Q A is close to the uncoupled horizontal tune Q x , and Q D is close to the uncoupled vertical tune Q y .The eigentunes can be determined with good precision from turn-by-turn data by filtering, Fourier transformation, and interpolation [13].For easier notation we will also use the quantities A;D 2Q A;D , and the symplectic conjugate R ÿSR T S of a matrix R, where S is the symplectic form and R T denotes the transposition of matrix R. For 2 2 and 4 4 matrices the symplectic forms are S 0 ÿ1 1 0 and S 0 B B @ 0 ÿ1 0 0 1 0 0 0 0 0 0 ÿ1 0 0 1 0 the complex conjugate is computed as We now note that a matrix G can be found that transforms the matrix M in Eq. ( 2) into an eigenbasis [11]: The underlined matrices denote one-turn maps in the eigenbasis.The matrix G and its inverse G ÿ1 can be expressed as [11] G gI ÿR R gI and G ÿ1 gI R ÿR gI ; (7) where I denotes the identity matrix.The quantities g and R can be computed as [11] g and The minimum tune approach is [11] Q min detjC Bj p sin A sin D (10) and global decoupling amounts to manipulations that result in detjC Bj 0. Note that the sign of detjC Bj is negative on sum resonances and positive on difference resonances [11].In Eqs. ( 8) and ( 9) the denominators become small when Q A ÿ D =2 !0, i.e., the globally decoupled state is approached.Under these circumstances we get with Eq. ( 10) since Q min Q.Thus g in Eq. ( 8) approaches 1.By inspection of Eqs. ( 6) and (7) it can be seen that lim Q!0 R 0. For weak coupling, the elements of the matrices B and C can be more than an order of magnitude smaller than those of the matrices A and D. As an example, consider a linear one-turn map of a RHIC model.In this model a single skew quadrupole of strength k 0:0022 m ÿ1 was inserted in an otherwise uncoupled linear machine, at a location with lattice functions x 11:8 m and y 44:9 m. k denotes the inverse focal length.The resulting minimum tune split is [1,10] and the one-turn map at the beginning of an arc, different from the location of the skew quadrupole, is ÿ1:87 48:88 ÿ0:07 ÿ0:49 ÿ0:11 2:33 0:00 ÿ0:04 0:01 ÿ0:54 0:80 12:90 0:00 ÿ0:02 ÿ0:11 ÿ0:52 The coupling information is in the small values of the offdiagonal submatrices which often sink into the noise floor when real data are used to fit a one-turn map.
With coupled motion energy is exchanged between the transverse planes with the beat frequency.To observe the maximum effect of the energy transfer from one plane to the other, one has to wait for half the beat period, or N turns with This can be seen in Figs.2(a) and 2(c).If we now consider a N-turn map and choose the turn number N with Eq. ( 14), the signalto-noise ratio of the coupling information is maximized.
The coupling information can also be obtained from the N-turn map.From Eq. ( 6) we get and from Eq. ( 9) which leads to Q min can be computed with Eq. ( 10).The signal-tonoise ratio of the matrix elements in C N B N is maximized when the factor is minimized.Often one can improve the results by minimizing this quantity through small variations in N from the start value given by Eq. ( 14).Based on the obtained value for C B, skew correctors can be set so that detjC Bj 0 when the correctors are included.This will be shown in Sec.V.

III. CONSTRUCTION OF TURN-BY-TURN 4-VECTORS
Assume that the beam position was measured at n monitors in a contiguous region over a number of turns.If the region contains only very small coupling sources, the horizontal and vertical orbits are decoupled for the passage of the beam during every turn.Such a region is typically an arc of a storage ring (see Fig. 1).The position and slope at the beginning of the arc in either plane can then be fitted from the beam positions at the n monitors.In principle two monitors are enough to obtain a fit result.However, beam position monitors may fail and by using more than two monitors the fit can be made resistant against a limited number of such failures.Monitors may deliver a status flag indicating whether the reported data are useful (see Refs. [14,15] for such an implementation at RHIC).Software checks can be implemented to reject unphysically large values (for example, see Ref. [16]).Furthermore, by using more monitors, the effect of random errors in the measured beam position on the fit result is reduced through averaging.
We denote by z i ; z 0 i T the beam position and slope in either plane at monitor i, relative to the closed orbit.The subscript 0 is used for the observation point at the beginning of the arc.The linear transfer map from the observation point to monitor i is [17 where i ; i denote the lattice functions and i the phase advance between the observation point and monitor i. Lattice functions and phase advances can be measured but in many cases values from the ideal lattice will be sufficient.If the beam has position z 0 and slope z 0 0 at the observation point, the position at monitor i is By minimizing the function the variables z 0 and z 0 0 can be fitted.To find the minimum of 2 z 0 ; z 0 0 we set the partial derivatives with respect to z 0 and z 0 0 to zero: Introducing the quantities In practice the coefficients M i 11 and M i 12 can be computed with Eq. ( 21) from model data, and the sums (26) can be calculated and stored in memory.When new data arrive, the closed orbit is subtracted from the measured positions and the sums (25) are recomputed while some beam position monitor data may be disregarded.The solution z 0 ; z 0 0 T is obtained with Eqs.(28).The procedure can be performed for the horizontal and vertical planes individually and the two 2-vectors are combined into a 4-vector x 0 ; x 0 0 ; y 0 ; y 0 0 T .This can be done for every turn that is recorded and thus a turn-by-turn 4-vector at the observation point is constructed.
The algorithm has only one division, by S in Eqs.(28).The value of S can be checked at the beginning.However, in a regular arc ( i arc const; i arc const) S is only zero if the phase advance between any two of the n monitors is .

IV. CONSTRUCTION OF A N-TURN MAP
We assume that m consecutive turns of a 4-vector z z x; x 0 ; y; y 0 T were fitted from turn-by-turn orbit data after a transverse kick.For the N-turn map one has To fit the matrix elements of M N the function is minimized.Note that in Eq. (30) weights w ij could be used for each addend corresponding to a fit parameter M ij .We have, however, chosen N such that the matrix elements of M N are not different by orders of magnitude and can avoid this complication.To find the minimum of 2 M N the partial derivatives with respect to the fit parameters M N ij are set to zero: We now introduce the two 4 4 matrices S a and S b with elements with which the condition (31) can be written as The direct solution for M N is M N S a S b ÿ1 : (34) Obviously, a solution for M N exists only if S b has an inverse.To avoid problems in an implementation, the condition detjS b j Þ 0 can be checked to ensure that a solution exists.The direct solution of Eq. (33) may not be the best way to solve the problem numerically.Better algorithms include LU decomposition with backsubstitution, Cholesky decomposition, and Gauss-Jordan elimination [18].The variances of the matrix elements M N ij can be obtained as [18] 2 M N ij The matrix elements M N ij are From this the variances are obtained as For hadron machines the matrix M N must be symplectic but the fit result (34) will be, in general, only nearly symplectic due to measurement errors in the beam position monitors and the effects of filamentation.Thus by symplectifying the matrix obtained with Eq. ( 34) the accuracy of the fit is likely to be improved.This can be done in the following way [19,20].A symplectic matrix M can be written as M expSU where U is a symmetric matrix.The matrix M can then be rewritten as where I is the identity matrix, and W is symmetric if and only if M is symplectic.Given a nearly symplectic matrix M ns , a matrix can be defined from which the exactly symmetric matrix W V V T =2 is obtained.The symplectified matrix M is then calculated using Eq. ( 40).Efficient symplectification routines are available [21].

V. MEASUREMENT AND CORRECTION OF LINEAR COUPLING
We assume that the eigentunes Q A;D were obtained from a Fourier transform of turn-by-turn data, and the matrix from a N-turn map.For a coupling measurement, the minimum tune approach Q min can then be calculated using Eq. ( 10).With Eq. ( 39) statistical errors M N ij for the matrix elements M N ij can be estimated, and with Eq. ( 19) the errors C ij and B ij for the matrix elements C ij and B ij can be obtained, respectively.An estimate for the error in the minimum tune approach can be computed from these quantities through error propagation as Large computed errors Q min indicate a significant mismatch between observed turn-by-turn orbit data and expectations from the model.This can happen with failures in the orbit acquisition system (for example, timing problems or wrong gains).Large errors can be used to reject coupling measurement results.
For a correction algorithm we work in a coordinate system, in which the linear motion is represented by circles in phase space.The transformation into the new coordinate system is provided by the matrix (46) and similar for B y [11].The matrix B is computed at the observation point.The matrices C and B in the old and new coordinate system are related via and the matrix K K can be written as We denote by i x the horizontal phase advance from the observation point to the skew quadrupole i, and use in the following i i x x ÿ i x (49) and In the new coordinate system we have for a number of weak skew quadrupoles [11] C C and The strength k i is the inverse focal length f i of skew quadrupole i. Assume a corrector family has the same skew corrector strength k 1 in all correctors and results in where the elements in C C 1 and B B 1 were divided by k 1 .For a global coupling correction we want to minimize the quantity From it follows with In principal two correctors or families are sufficient to correct linear coupling globally (unless their matrices K K 1 and K K 2 are linear dependent).After the first strength k 1 has been found, the second strength k 2 can be found with Eqs. ( 56) and (57), by replacing If there are more than two corrector families, their effectiveness e i can be tested by computing where Q min;i the minimum tune approach for K K k i K K i .After the most effective corrector k m is set, the test can be repeated for the next corrector with K K replaced by K K k m K K m , and so on.The corrector strengths can also be distributed among the available correctors using other constraints.

VI. APPLICATION AT THE RELATIVISTIC HEAVY ION COLLIDER
The above described algorithm for the global linear coupling correction has been implemented within the RHIC injection optimization application.After the beam is injected, turn-by-turn data are automatically acquired from 12 beam position monitors in the horizontal and 12 monitors in the vertical plane.1024 turns are recorded in each of the beam position monitors while the beam exhibits injection oscillations.The application has access to an on-line machine optics model and can read and set skew corrector strengths [22].
The implementation was first tested with simulated turn-by-turn data and defined skew quadrupole errors in the optics model.Then, using turn-by-turn data from machine operation, the Q min computed from the N-turn maps was compared with the Q min obtained by bringing the tunes together.This is shown in Fig. 3.In part (a) the measured tunes are plotted as a function of the horizontal tune set point.From this data one obtains Q min 0:0137.In Fig. 3(b) the Q min computed from N-turn maps is shown for the same horizontal tune set points, along with the errors computed with Eqs.(39) and (43).Over the scanned horizontal tune set points, the Q min from the N-turn maps differs from the Q min obtained by bringing the tunes together by less than 0.002 in all but two cases.The coupling measurement returns useful values over a range of tune differences jQ x ÿ Q y j sufficient to cover operational needs (both RHIC tunes are kept between 0.2 and 0.25).The computed errors appear to be a conservative estimate.They become large for large tune differences.
As an example for a coupling measurement and correction we show the first test of the algorithm in operation, performed in the RHIC Blue ring with a deuteron beam.In Figs.2(a Here too the coupling is visible with two large peaks in both spectra.From the measured eigentunes, the optimum turn number N is determined with Eq. ( 14) and the subsequent minimization of the quantity (20) through small variations in N. We get N 67 and the fitted N-turn map is The result of the coupling correction is shown in Figs.2(e) -2(h).In 2(e) and 2(g) the reduction in beating due to coupling is recognizable.The recoherence after 650 turns is due to synchrotron motion and nonzero chromaticity [23].650 turns are a synchrotron period.The reduction in coupling is also visible in the spectra in Figs.The Q min reached by the correction was confirmed by a tune measurement after moving the tunes together.
In operation it was found that the coupling correction cannot improve Q min significantly beyond 0.002.This is consistent with the Q min 0:0011 predicted for the next correction.Analyzing about 100 injections, it was found that detjC Bj was negative in 10% of the cases, typical of sum resonances [11].These cases are rejected for correction.With increasing chromaticities the signal decoheres faster.Although a correction can still be computed with a few hundred turns, the part during which the signal is decohered has to be filtered out.It was also found that the computed errors for Q min are smaller for small betatron oscillations.This indicates that nonlinearities at large amplitudes, not included in the model, lead to discrepancies between the measured turn-by-turn data and expectations from the model.

VII. SUMMARY
A method for global coupling measurement and correction is presented that is based on N-turn maps fitted from turn-by-turn beam position data.N is chosen so as to maximize the signal-to-noise ratio of the coupling information.By using more than two monitors per plane the robustness is increased and the effect of random errors in beam position monitors is ameliorated.No tune change is needed for either the measurement or the correction.The method is implemented for operation as part of the injection optimization application at the Brookhaven Relativistic Heavy Ion Collider.It allows a global coupling correction within seconds after the turn-by-turn data are acquired.The method may be used in situations other than injection, when turn-by-turn data of free beam oscillations can be acquired.

FIG. 1 .
FIG.1.Observation point at the beginning of an arc and beam trajectory in the arc.

PRST 3 FIG. 2
FIG. 2. (Color) Turn-by-turn signals and spectra before and after a coupling correction.In (a) and (c) the horizontal and vertical injection oscillations are shown before a coupling correction.In (b) and (d) the respective spectra are depicted.(b) also shows the computed Q min , its error, and the predicted Q min after a coupling correction.(e) -(h) show the situation after the computed corrector values were applied.
FIG. 3. (Color) Measured tunes as a function of the horizontal tune set point in (a).(b) shows the Q min obtained from N-turn maps, and the computed errors, for the same tune set points.Also shown in (b) is the Q min that can be extracted from the tune measurement depicted in (a).
which Q min 0:0064 is obtained.The error according to formula (43) is 0.0028.RHIC has three families for global decoupling, due to the sixfold symmetry of the machine.With the algorithm, two of those families are selected to minimize the coupling.The predicted Q min after correction is 0.0003 [see Fig.2(b)].A nonzero prediction is a sign of measurement errors in the N-turn map, or a mismatch between the optics model and the machine.The correction can be implemented by pressing a single button in the application.
2(f) and 2(h).After correction, we have is reduced to 0.0023 with an error of 0.0006.Note the reduction in the matrix elements M 14 and M 32 .