- Email: [email protected]
Fast retrace of particle orbits through a large, inhomogeneous magnetic field
NUCLEAR
INSTRUMENTS
AND
METHODS
II 4
0974) 3 2 1 - 3 2 6 ;
NORTH-HOLLAND
©
PUBLISHING
CO.
FAST RETRACE OF PARTICLE ORBITS T H R O U G H A LARGE,
I N H O M O G E N E O U S MAGNETIC FIELD McKEE
R.J.
Division of" Physics, National Research Council of Canada, Ottawa, Ontario, K1A OR6, Canada Received 8 June 1973 A fast m e t h o d of retracing orbits t h r o u g h a large, i n h o m o g e n e o u s magnetic field has been developed for applications in particlephysics experiments. The m o m e n t a of particles are found by simple linear interpolation in a mesh o f previously calculated representative orbits. T h e essential feature o f the scheme is to perform a t r a n s f o r m a t i o n on the set of observables defining the orbit before doing the interpolation. This is done to reduce both
interpolation errors a n d the required n u m b e r o f mesh points. A detailed development of the m e t h o d is presented for an actual apparatus. For this a p p a r a t u s particles were retraced in 0.6 m s / orbit using a C D C 6600 computer, a n d for 100 mesh points interpolation errors were found to be smaller t h a n experimental errors by more t h a n an order o f magnitude.
1. Introduction
in ref. 2. In the next section the retrace mesh formalism is developed for our apparatus.
During the course of a particle-physics experiment 1) that uses a spectrometer consisting of a bending magnet and multiwire proportional chambers, we have developed a method of rapid and accurate orbit retracing through the large, inhomogeneous magnetic field. The method is to derive a relatively simple, approximate function relating the momentum of the particle to the observables defining the orbit (in our case wire chamber coordinfites). The function is then used to find the momenta of experimental data. Lechanoine et al. 2) have also developed such a method. They compute the m o m e n t u m for a few representative orbits and use the results to define a fimction that gives the momentum explicitly in terms of the observables; the functional form they develop consists of sums of products of Tchebycheff polynomials fitted to the calculated orbits. In our method we first of all define a new set of variables in terms of the observables; using the momenta calculated from a few representative orbits, we develop a function that gives the momentum in terms of the new variables instead of directly in terms of the observables. By choosing the representative orbits so that they form a rectangular mesh in the space defined by the new variables, the procedure we use is linear interpolation using the momenta at the nearest mesh points. The whole point of using a new set of independent variables instead of the observables is to allow us to easily construct a rectangular mesh and to keep the interpolation errors small f o r a small number of mesh points. It turns out that the use of our method to analyze data is ten times faster than the method quoted 321
2. The retrace mesh
Fig. 1. shows a simplified drawing of the apparatus. It consists of a target, bending magnet, and two spectrometer arms, Each arm consists of four wire planes; two X views (horizontal) and two Y views (vertical). The X2 and Y2 wire planes are located Y *,
ARM
ONE
i
i
/
_/
W C = X - AND Y VIEW WIRE CHAMBERS
-~',
I :'-
MAGNET TIPSi
POLE
l/ WC I
~
" ~ ~ ~//~
X
TOP VIEW
WC 2
. . . . .
PLANE
~L . . . .
WC L
SIDE ALONG
WC 2
VIEW A WIRE
CHAMBER
ARM
Fig. 1. A simplified drawing of the a p p a r a t u s for which the retrace m e s h w a s developed.
322
R.J. McKEE
downstream from the X1 and Y~ planes. The spectrometer is specifically designed to study reactions involving just two final charged particles whose orbits are indicated by the dashed lines. However, in what is to follow, each arm is treated separately and in the end has its own retrace mesh. Since the four wire coordinates of an arm completely specify the orbit downstream of the magnet, we seek a function of four independent variables relating the momentum at the target to the four wire coordinates. But first it is necessary to define exactly what constitutes a calculated orbit retrace for this apparatus; with no positional measuring counters between target and magnet (the small, point-like target is actually quite deep in the fringing field), the target itself must be used to complete the specification of an orbit. We assume that the horizontal projection of the orbit passes through the center of the target; for any arbitrary but reasonable set of four wire chamber coordinates there is no guarantee that the orbit itself passes through the center of the target. We then guess a momentum magnitude and trace the orbit back through the magnet by numerical integration of the equation of motion using the previously obtaiaed map of the magnetic field. If the orbit projection misses the center of the target, the momentum magnitude is systematically varied, by Newton's method, until the orbit meets our criteria. We therefore want a simple function of four inde-
pendent variables which approximates the above calculation. We have adapted the following: A 4-dimensional rectangular mesh of representative orbits is chosen, and the momentum components at the target for each mesh point are calculated by the long procedure described above, Each mesh point is in a 4-dimensional space since each orbit is defined by the four wire coordinates. The momentum components of experimental data are found by linear interpolation in the 4-dimensional mesh. The problem is the selection of the mesh points. We wish to keep the total number of mesh points small, but the points must be close enough together so that interpolation errors are also kept small. [t would be difficult to make the mesh directly in terms of the wire coordinates. This can be seen by studying fig. 2. The diagram shows a typical scatterplot in one of the arms of wire coordinates X z vs Xj, where X 2 is the horizontal view downstream from X1. The shaded area contains the wanted experimental data. With the strong correlation between Xj and X 2 it is far better to construct the mesh in terms of the variables V1 and V2 shown as grid lines on the diagram with the points of intersection as the mesh points. The problem then reduces to defining a suitable mapping between the X's and the V's. The problem can be stated more generally. We seek a set of new variables
V k = F k ( X 1 , X 2 , Y1, Y2),
Vr
A
Xi g(X2)
X~ Fig. 2. A plot of wire coordinates X2 vs X1 for one spectrometer arms. The wanted data is represented shaded region. The new variables V1 and V2 are shown lines superimposed over the region of correlation. For the curvature of g(Xz) is exaggerated.
of the by the as grid clarity,
k = l ..... 4,
(1)
where the transformation functions Fk are simple, have a unique inverse, reduce interpolation errors, and for a rectangular region in V-space, cover the region of orbits allowed by the spectrometer arm. Since the basic idea is to calculate the momentum directly as a function of the V's instead of the observables, a unique inverse of eq. (1) is needed, at least over the range of the mesh; that is, we start with a mesh point in V-space and then we must be able to calculate the four wire coordinates uniquely so that the momentum can be found by numerical integration of the orbit. Since the V's themselves have no physical significance, they can be scaled in any convenient manner. After choosing the mesh ranges along each Vk axis, the transformation functions Fk are then chosen so that the allowed orbits nearly fill the rectangular volume thus defined. The Fk are also chosen so that the momentum components are nearly linear in the V's; this reduces interpolation errors, or, equivalently, reduces the number of mesh points needed to obtain a given precision. The transformations are also constructed to be relatively simple
FAST RETRACE
OF P A R T I C L E
since each experimental event must be transformed by eq. (1) before interpolation. 'With the help of fig. 2 these ideas are now applied to derive the transformations for V1 and V2. For simplicity, we define these two transformations to be independent of Y~ and Y2. In order to reduce the interpolation error in finding, say, the x-component of momentum Px (see fig. 1 for the definition of the coordinate axis), it is desirable to define V~ in terms of X~ and 2(2 such that it is linear, or nearly so, in px. This property is represented on the diagram by the uneven spacing of the vertical lines of constant V~. For simplicity if we require linearity in Px for only one value of V2, then Va is a function of just one of the horizontal wire coordinates. This comes about as fo]tlows: Consider the orbit restricted to the median plane of the magnet (fig. 1). As the lab angle of the reaction is varied and the orbit traced forward from the center of the target to the wire chambers, a one-toone correspondence between px, X~, and )(2 is generated; treating X2 as the independent variable we get two functions: Px = f ( X 2 ) ,
(2)
x, = g(Xg.
(3)
The second relationship, eq. (3), is represented in fig;. 2 by the central line; reactions outside the median plane trace other lines, more or less parallel, through the shaded region. We simply define V1 proportional to f(X2). Since the central line on the diagram is, by definition, a line of constant V2, we define V2 proportional to X 1 - g ( X 2 ) . The mapping function for V1 and Vz, given essentially by f(X2) and g(X2), can only be found numerically. However, it is adequate to approximate each function by a quadratic polynomial. This is done by calculating p~, X~, and X2 for three different lab angles. The final form of the transformation is V1 = ao+aLX2+a2X
(4)
2,
V2 = X , - ( b o + b, X 2 + b2 XZ2).
(5)
If we let the superscript I stand for the three different lab angles, we find the coefficients of eq. (~,) by fitting a quadratic to V(1)
X(3) 2
--A2
t'x
-- ~x
_(l)
_(1)X
X~2n vs X~ ) + n~3) #1) (P~ - v ~ J,
(6)
and the coefficients of eq. (5) by fitting another quadratic to X~2n vs X[ n. With X~21)and X2(3) chosen to be slightly beyond the limits of the spectrometer arm
323
ORBITS
acceptance (with X~22) near the center), there is no reason to extend the mesh further. By eq. (6), V~ = X2 at these two end points. Hence, the form of eq. (6) automatically gives us the mesh range along the VI axis. Because of the strong correlation between X1 and X 2, the mesh range along the V2 axis is small and can be chosen somewhat larger than the thickness of the shaded region in the X~ direction (fig. 2). There are some problems that may arise, however. As noted earlier, a unique inverse transformation must exist over the range of the mesh. Solving for X 2 in eq. (4) results in X2 = a, [__(l+q)~
1] '
r/C0,
2a 2 l/"1 - - a
O
, = 0,
(7)
a!
where
q
--
4az(V1 - a 0 ) a2
(8)
If f (X2) has an extremum over the mesh range, then both signs are possible and no unique inverse exists. In this case we define V1 = XE. With no extremum the plus sign is used if V1 monotonically increases with X2 and the minus sign for the decreasing case. Another problem occurs in the Jacobian peak region, if one exists. The function g(X2) sharply doubles back on itself forming a hook at one end. In that case we choose the larger of the two branches for the mesh range in V1 and increase the mesh range in V2 to cover the hook region. While the transformation just described for V1 is suitable for Px by construction, large interpolation errors may result if it is used for the mesh in py, the y-component of momentum. A different transformation is used for py; its definition is analogous to the one already described for Px. All that need be done is to replace p~) by if}) in eq. (6) thus giving a different set of coefficients for eq. (4). In our apparatus the vertical component of momentum p~ is always small compared to the horizontal component. For simplicity, no new transformation is used for Pz. Since the momentum in the center of momentum system (c.m.s.) varies less rapidly than in the lab system, our meshes are actually calculated for momentum components in the c.m.s. This also speeds up final data analysis since most of the quantities of interest are also in the c.m.s. We now consider the transformations of the vertical wire coordiantes. A scatterplot of Y2 vs Y1 looks quite
324
R.J. McKEE
similar to fig. 2. Hence, the transformations of Y, and Y2 to V3 and V4 are analogous to the two transformations of the horizontal wire coordinates. But for our apparatus there are two simplifications. The geometry restricts the orbits close to the median plane. Thus, Px and py are very insensitive to the vertical wire coordinate ; Pz varies quite linearly with Yz already and is small. So we define V3 proportional to Y2 immediately without using a fitted polynomial. The second simplification is due to the high degree of mirror symmetry in the magnetic field; as a result, the orbits themselves have mirror symmetry. This means that we are able to restrict the range of the mesh to include only orbits above the median plane. With the above remarks in mind, we now derive explicit transformations for V3 and V4 analogous to Vl and V2. Because of the mirror symmetry it is convenient to define vertical wire coordinates relative to the median plane: Y( = Y I - Ym,
Y2 = ]/2- Ym,
(91
where Ym, the wire coordinate in the median plane, is assumed to be the same for both wire planes. First of all, no function analogous to f(X2) is needed. The function analogous to g(X2) is found by restricting the lab angle to a fixed value (near the center of the setting for the spectrometer arm) and varying the plane of the reaction. As we trace the orbit forward, we get a functional relationship between Y, and Y2 which we write as
,20
Y; = h(Y~).
(10)
Analogous to the definition of V 2 , w e define V4 proportional to Y [ - h ( Y ~ ) . But by mirror symmetry h(Y~) is an odd function' both Y~ and Y~ change signs simultaneously. Therefore, the polynomial fitted to h(Y~) must also be odd. The final forms for the last two transformations are ~3 = Y2, I~k = Y ; - Y [ ( C o + C ,
(11) Y ~ + C 2 Y~").
(12)
The coefficients are found by calculating Y1 and Y2 for three different (non-zero) values of the angle of the reaction plane and fitting y:~2vs Y[/Y~ to a quadratic. As mentioned before, the mesh range along the V3 axis is restricted to positive values. The mesh range for V4 is small and corresponds to the thickness of the Y2 v s YI correlation in the direction of Yr. Also note that if Y~<0 for an experimental data event, the signs of both Y~ and Y~ must be changed before applying eqs. (11) and (12). The resulting sign found forp~ after interpolation must then be changed before any further analysis is done for the event. In addition to the meshes in the three m o m e n t u m components, we have found it useful to include a fourth mesh in z, the vertical distance by which the orbit itself misses the center of the target. This quantity can be used to cut events not originating from the target. The previous remarks made for p: also apply for z. Specifically, no new transformations are needed for z; also, if Y ; < 0 , the sign of z found after interpelation must be changed.
Ioo
z
ARM TWO p~55 GeV/c
ARM ONE p~17 GeV/c
80 [
w
~ ~o cr Z
40
z
© ~-O 15
,
-OIO
-005
O
005
OtO
Oi5 -oi5
-oio
-005
O
O 05
O '10
rl OI5
I00 (p- p ' ) / p
Fig. 3. Histogram of interpolation error (in %) for 1000 Monte-Carlo events, p = original lab momentum, p ' = interpolated momentum.
FAST RETRACE
OF P A R T I C L E
To conclude this section we present the linear interpolation equation, which is formulated to minimize the number of arithmetic operations. In this formulation the values of the momentum components and z at each mesh point are not used directly; instead certain linear combinations of the values at 16 neighboring mesh points are used and must be calculated in aclvance. To simplify things, let P stand for one of the three momentum components or z. After transforming an experimental event given by the four wire coordinates, we must find into which 4-dimensional rectangle the data point falls. Let Vk0 and A V k be the first mesh point and the (constant) mesh interval along the V k axis. Let I k --
v~-V~o AV~
Xk = t k - - [ t k ] ,
(13)
where [t] is the largest integer less than or equal to t. T]he four [tk] show which rectangle the data point is in while the Xk give the relative position inside the rectangle. We are assuming that the data point actually falls somewhere inside the pre-established mesh. Let P~B~6be the value of P at one of the 16 corners of the 4-dimensional rectangle; the kth subscript value is zero or one depending on whether P~a~o is on the low or high side of the rectangle along the Vk axis. The 4-dimensional linear interpolation is then given by P = C1 + C2 xl
+ ( C 3 + C 4 X l ) x 2 -~-
+ [C5 + C6 x~ +(C7 + Cs xl)x2] x3 + -'F'-{ C 9 + C l o X 1 -~- ( C l i -{- C 12 XI)X2 -~-
-+-[C13-[-C14x I + ( C I 5 + C I 6 X l ) X 2 3 X 3 } X 4 ,
(14)
where the Ci are given in table 1.
3. Speed and accuracy A typical application uses a mesh that is 5 by 5 by 2 by 2 for V~ through I/4, which amounts to only 100 mesh points and sixteen 4-dimensional rectangles. The small number of mesh points is an important consideration since much computer time is needed to calculate the momentum for each mesh point (about ¼ of a second on the CDC 7600). Also note that the space requirements for tLe mesh in the computer is small; with 16 rectangles, 1024 memory locations are needed for each spectrometer arm. The computing time needed to process data for each spectrometer arm is 0.6 ms/event (on the C D C 6600). This includes the time needed to transform the four
325
ORBITS TABLE 1
The 16 i nt e rpol a t i on c ons t a nt s used in eq. (14).
i 1 2 3 4 5 6 7 8 9 lO Il 12 13 14 15 16
G Poooo Plooo -- C1 Poloo - C1 P l l o O - (C1 + C 2 + C3) Poolo - C1 PlOlO -- (C1 + C2 + C5) P011o - (C1 + Ca + C5) 7 Pl110 - ~ C, l=l PO001 -- C1 PlOOl-(CI+C2+C9) POlOl - (C1 + C3 + Co)
Pl101-(C1+C2+C3+C4+C9+C10+Cll) POOll - (C1 + C5 + C9) P l O l l - (C1 + C2 + C5 + C6 + C9 + Clo + ClZ) Pol 11 - (C1 + C3 + C5 + C7 + C9 + C11 + C13) 15 P l l l l -- t__~lC~
wire coordinates and to interpolate all three momentum components and z. Lechanoine et al. 2) quote times of 6 and 15 ms/event (also on a CDC 6600) for two different cases that they describe. The difference in speed is not hard to find. Their use of a function of five observables consisting of sums of products of Tchebycheff polynomials is really a more complicated interpolation scheme in a five dimensional mesh of calculated orbits. (The mesh points are chosen according to certain properties of Tchebycheff polynomials.) Hence, a large number of mesh points are used in their calculation while we use only the nearest 16 points. With 100 mesh points the error introduced by the method is less than 0.5% which compares favorably with the accuracy quoted in ref. 2. A Monte-Carlo simulation was done to check the accuracy. Events chosen uniformly over a typical spectrometer setting were traced forward by numerical integration and retraced by the mesh. A histogram of the interpolation errors in both arms is shown in fig. 3. The error is over an order of magnitude smaller than the experimental resolution of the spectrometer. The biggest source of error occurs in the interpolation along the V z axis due to the rapid variation of momentum with Va for a fixed VI. This is why we generally use as many mesh peints along V 2 as along V1. Although never tried, the accuracy may be improved by rotating the lines of constant V~ so that they end up perpendicular to the central line X~ = g ( X E ) (fig. 2). Even though in
326
R.J.
our example the n u m b e r o f mesh points along V3 and V4 are as few as possible, i n t e r p o l a t i o n errors along these two axes are the smallest o f all; this is due to the already mentioned fact that the horizontal m o m e n t u m c o m p o n e n t is insensitive to V 3 and //4 while Pz and z are nearly linear in the vertical wire coordinates. 4. Conclusion To summarize: A multi-dimensional rectangular mesh is chosen in a space o f new variables whose dimensionality is equal to the n u m b e r of observables needed to define an orbit. A suitable set o f transformations between the new variables and the observables is chosen so that the mesh covers the physically meaningful orbits and linearizes at least one momentum component. Other t r a n s f o r m a t i o n s may be used to linearize other m o m e n t u m components. Once the m o m e n t u m c o m p o n e n t s are calculated at each mesh point prior to d a t a processing, the m o m e n t u m of
McKEE
each experimental event is found by t r a n s f o r m i n g the observables into the new variables and linearily interpolating in the multi-dimensional mesh. The detailed d e v e l o p m e n t of section 2 was presented to show how one might go a b o u t applying the idea to a specific apparatus. The description in section 2 might suggest that this retrace m e t h o d is strongly tied to a specific apparatus. But the basic idea - suitable t r a n s f o r m a t i o n o f observables and linear interpolation in a mesh f o r m e d by the t r a n s f o r m e d variables to find the m o m e n t u m quickly and accurately - certainly has general application.
References ~) H. L. Anderson, D. Larson, L. Myrianthopoulos, L. Dubal, C. K. Hargrove, E. P. Hincks, R. J. McKee, H. Mes, D. Kessler, and A. C. Thompson, to be published. 2) C. Lechanoine, M. Martin, and H. Wind, Nucl. Instr. and Meth. 69 (1969) 122.
INSTRUMENTS
AND
METHODS
II 4
0974) 3 2 1 - 3 2 6 ;
NORTH-HOLLAND
©
PUBLISHING
CO.
FAST RETRACE OF PARTICLE ORBITS T H R O U G H A LARGE,
I N H O M O G E N E O U S MAGNETIC FIELD McKEE
R.J.
Division of" Physics, National Research Council of Canada, Ottawa, Ontario, K1A OR6, Canada Received 8 June 1973 A fast m e t h o d of retracing orbits t h r o u g h a large, i n h o m o g e n e o u s magnetic field has been developed for applications in particlephysics experiments. The m o m e n t a of particles are found by simple linear interpolation in a mesh o f previously calculated representative orbits. T h e essential feature o f the scheme is to perform a t r a n s f o r m a t i o n on the set of observables defining the orbit before doing the interpolation. This is done to reduce both
interpolation errors a n d the required n u m b e r o f mesh points. A detailed development of the m e t h o d is presented for an actual apparatus. For this a p p a r a t u s particles were retraced in 0.6 m s / orbit using a C D C 6600 computer, a n d for 100 mesh points interpolation errors were found to be smaller t h a n experimental errors by more t h a n an order o f magnitude.
1. Introduction
in ref. 2. In the next section the retrace mesh formalism is developed for our apparatus.
During the course of a particle-physics experiment 1) that uses a spectrometer consisting of a bending magnet and multiwire proportional chambers, we have developed a method of rapid and accurate orbit retracing through the large, inhomogeneous magnetic field. The method is to derive a relatively simple, approximate function relating the momentum of the particle to the observables defining the orbit (in our case wire chamber coordinfites). The function is then used to find the momenta of experimental data. Lechanoine et al. 2) have also developed such a method. They compute the m o m e n t u m for a few representative orbits and use the results to define a fimction that gives the momentum explicitly in terms of the observables; the functional form they develop consists of sums of products of Tchebycheff polynomials fitted to the calculated orbits. In our method we first of all define a new set of variables in terms of the observables; using the momenta calculated from a few representative orbits, we develop a function that gives the momentum in terms of the new variables instead of directly in terms of the observables. By choosing the representative orbits so that they form a rectangular mesh in the space defined by the new variables, the procedure we use is linear interpolation using the momenta at the nearest mesh points. The whole point of using a new set of independent variables instead of the observables is to allow us to easily construct a rectangular mesh and to keep the interpolation errors small f o r a small number of mesh points. It turns out that the use of our method to analyze data is ten times faster than the method quoted 321
2. The retrace mesh
Fig. 1. shows a simplified drawing of the apparatus. It consists of a target, bending magnet, and two spectrometer arms, Each arm consists of four wire planes; two X views (horizontal) and two Y views (vertical). The X2 and Y2 wire planes are located Y *,
ARM
ONE
i
i
/
_/
W C = X - AND Y VIEW WIRE CHAMBERS
-~',
I :'-
MAGNET TIPSi
POLE
l/ WC I
~
" ~ ~ ~//~
X
TOP VIEW
WC 2
. . . . .
PLANE
~L . . . .
WC L
SIDE ALONG
WC 2
VIEW A WIRE
CHAMBER
ARM
Fig. 1. A simplified drawing of the a p p a r a t u s for which the retrace m e s h w a s developed.
322
R.J. McKEE
downstream from the X1 and Y~ planes. The spectrometer is specifically designed to study reactions involving just two final charged particles whose orbits are indicated by the dashed lines. However, in what is to follow, each arm is treated separately and in the end has its own retrace mesh. Since the four wire coordinates of an arm completely specify the orbit downstream of the magnet, we seek a function of four independent variables relating the momentum at the target to the four wire coordinates. But first it is necessary to define exactly what constitutes a calculated orbit retrace for this apparatus; with no positional measuring counters between target and magnet (the small, point-like target is actually quite deep in the fringing field), the target itself must be used to complete the specification of an orbit. We assume that the horizontal projection of the orbit passes through the center of the target; for any arbitrary but reasonable set of four wire chamber coordinates there is no guarantee that the orbit itself passes through the center of the target. We then guess a momentum magnitude and trace the orbit back through the magnet by numerical integration of the equation of motion using the previously obtaiaed map of the magnetic field. If the orbit projection misses the center of the target, the momentum magnitude is systematically varied, by Newton's method, until the orbit meets our criteria. We therefore want a simple function of four inde-
pendent variables which approximates the above calculation. We have adapted the following: A 4-dimensional rectangular mesh of representative orbits is chosen, and the momentum components at the target for each mesh point are calculated by the long procedure described above, Each mesh point is in a 4-dimensional space since each orbit is defined by the four wire coordinates. The momentum components of experimental data are found by linear interpolation in the 4-dimensional mesh. The problem is the selection of the mesh points. We wish to keep the total number of mesh points small, but the points must be close enough together so that interpolation errors are also kept small. [t would be difficult to make the mesh directly in terms of the wire coordinates. This can be seen by studying fig. 2. The diagram shows a typical scatterplot in one of the arms of wire coordinates X z vs Xj, where X 2 is the horizontal view downstream from X1. The shaded area contains the wanted experimental data. With the strong correlation between Xj and X 2 it is far better to construct the mesh in terms of the variables V1 and V2 shown as grid lines on the diagram with the points of intersection as the mesh points. The problem then reduces to defining a suitable mapping between the X's and the V's. The problem can be stated more generally. We seek a set of new variables
V k = F k ( X 1 , X 2 , Y1, Y2),
Vr
A
Xi g(X2)
X~ Fig. 2. A plot of wire coordinates X2 vs X1 for one spectrometer arms. The wanted data is represented shaded region. The new variables V1 and V2 are shown lines superimposed over the region of correlation. For the curvature of g(Xz) is exaggerated.
of the by the as grid clarity,
k = l ..... 4,
(1)
where the transformation functions Fk are simple, have a unique inverse, reduce interpolation errors, and for a rectangular region in V-space, cover the region of orbits allowed by the spectrometer arm. Since the basic idea is to calculate the momentum directly as a function of the V's instead of the observables, a unique inverse of eq. (1) is needed, at least over the range of the mesh; that is, we start with a mesh point in V-space and then we must be able to calculate the four wire coordinates uniquely so that the momentum can be found by numerical integration of the orbit. Since the V's themselves have no physical significance, they can be scaled in any convenient manner. After choosing the mesh ranges along each Vk axis, the transformation functions Fk are then chosen so that the allowed orbits nearly fill the rectangular volume thus defined. The Fk are also chosen so that the momentum components are nearly linear in the V's; this reduces interpolation errors, or, equivalently, reduces the number of mesh points needed to obtain a given precision. The transformations are also constructed to be relatively simple
FAST RETRACE
OF P A R T I C L E
since each experimental event must be transformed by eq. (1) before interpolation. 'With the help of fig. 2 these ideas are now applied to derive the transformations for V1 and V2. For simplicity, we define these two transformations to be independent of Y~ and Y2. In order to reduce the interpolation error in finding, say, the x-component of momentum Px (see fig. 1 for the definition of the coordinate axis), it is desirable to define V~ in terms of X~ and 2(2 such that it is linear, or nearly so, in px. This property is represented on the diagram by the uneven spacing of the vertical lines of constant V~. For simplicity if we require linearity in Px for only one value of V2, then Va is a function of just one of the horizontal wire coordinates. This comes about as fo]tlows: Consider the orbit restricted to the median plane of the magnet (fig. 1). As the lab angle of the reaction is varied and the orbit traced forward from the center of the target to the wire chambers, a one-toone correspondence between px, X~, and )(2 is generated; treating X2 as the independent variable we get two functions: Px = f ( X 2 ) ,
(2)
x, = g(Xg.
(3)
The second relationship, eq. (3), is represented in fig;. 2 by the central line; reactions outside the median plane trace other lines, more or less parallel, through the shaded region. We simply define V1 proportional to f(X2). Since the central line on the diagram is, by definition, a line of constant V2, we define V2 proportional to X 1 - g ( X 2 ) . The mapping function for V1 and Vz, given essentially by f(X2) and g(X2), can only be found numerically. However, it is adequate to approximate each function by a quadratic polynomial. This is done by calculating p~, X~, and X2 for three different lab angles. The final form of the transformation is V1 = ao+aLX2+a2X
(4)
2,
V2 = X , - ( b o + b, X 2 + b2 XZ2).
(5)
If we let the superscript I stand for the three different lab angles, we find the coefficients of eq. (~,) by fitting a quadratic to V(1)
X(3) 2
--A2
t'x
-- ~x
_(l)
_(1)X
X~2n vs X~ ) + n~3) #1) (P~ - v ~ J,
(6)
and the coefficients of eq. (5) by fitting another quadratic to X~2n vs X[ n. With X~21)and X2(3) chosen to be slightly beyond the limits of the spectrometer arm
323
ORBITS
acceptance (with X~22) near the center), there is no reason to extend the mesh further. By eq. (6), V~ = X2 at these two end points. Hence, the form of eq. (6) automatically gives us the mesh range along the VI axis. Because of the strong correlation between X1 and X 2, the mesh range along the V2 axis is small and can be chosen somewhat larger than the thickness of the shaded region in the X~ direction (fig. 2). There are some problems that may arise, however. As noted earlier, a unique inverse transformation must exist over the range of the mesh. Solving for X 2 in eq. (4) results in X2 = a, [__(l+q)~
1] '
r/C0,
2a 2 l/"1 - - a
O
, = 0,
(7)
a!
where
q
--
4az(V1 - a 0 ) a2
(8)
If f (X2) has an extremum over the mesh range, then both signs are possible and no unique inverse exists. In this case we define V1 = XE. With no extremum the plus sign is used if V1 monotonically increases with X2 and the minus sign for the decreasing case. Another problem occurs in the Jacobian peak region, if one exists. The function g(X2) sharply doubles back on itself forming a hook at one end. In that case we choose the larger of the two branches for the mesh range in V1 and increase the mesh range in V2 to cover the hook region. While the transformation just described for V1 is suitable for Px by construction, large interpolation errors may result if it is used for the mesh in py, the y-component of momentum. A different transformation is used for py; its definition is analogous to the one already described for Px. All that need be done is to replace p~) by if}) in eq. (6) thus giving a different set of coefficients for eq. (4). In our apparatus the vertical component of momentum p~ is always small compared to the horizontal component. For simplicity, no new transformation is used for Pz. Since the momentum in the center of momentum system (c.m.s.) varies less rapidly than in the lab system, our meshes are actually calculated for momentum components in the c.m.s. This also speeds up final data analysis since most of the quantities of interest are also in the c.m.s. We now consider the transformations of the vertical wire coordiantes. A scatterplot of Y2 vs Y1 looks quite
324
R.J. McKEE
similar to fig. 2. Hence, the transformations of Y, and Y2 to V3 and V4 are analogous to the two transformations of the horizontal wire coordinates. But for our apparatus there are two simplifications. The geometry restricts the orbits close to the median plane. Thus, Px and py are very insensitive to the vertical wire coordinate ; Pz varies quite linearly with Yz already and is small. So we define V3 proportional to Y2 immediately without using a fitted polynomial. The second simplification is due to the high degree of mirror symmetry in the magnetic field; as a result, the orbits themselves have mirror symmetry. This means that we are able to restrict the range of the mesh to include only orbits above the median plane. With the above remarks in mind, we now derive explicit transformations for V3 and V4 analogous to Vl and V2. Because of the mirror symmetry it is convenient to define vertical wire coordinates relative to the median plane: Y( = Y I - Ym,
Y2 = ]/2- Ym,
(91
where Ym, the wire coordinate in the median plane, is assumed to be the same for both wire planes. First of all, no function analogous to f(X2) is needed. The function analogous to g(X2) is found by restricting the lab angle to a fixed value (near the center of the setting for the spectrometer arm) and varying the plane of the reaction. As we trace the orbit forward, we get a functional relationship between Y, and Y2 which we write as
,20
Y; = h(Y~).
(10)
Analogous to the definition of V 2 , w e define V4 proportional to Y [ - h ( Y ~ ) . But by mirror symmetry h(Y~) is an odd function' both Y~ and Y~ change signs simultaneously. Therefore, the polynomial fitted to h(Y~) must also be odd. The final forms for the last two transformations are ~3 = Y2, I~k = Y ; - Y [ ( C o + C ,
(11) Y ~ + C 2 Y~").
(12)
The coefficients are found by calculating Y1 and Y2 for three different (non-zero) values of the angle of the reaction plane and fitting y:~2vs Y[/Y~ to a quadratic. As mentioned before, the mesh range along the V3 axis is restricted to positive values. The mesh range for V4 is small and corresponds to the thickness of the Y2 v s YI correlation in the direction of Yr. Also note that if Y~<0 for an experimental data event, the signs of both Y~ and Y~ must be changed before applying eqs. (11) and (12). The resulting sign found forp~ after interpolation must then be changed before any further analysis is done for the event. In addition to the meshes in the three m o m e n t u m components, we have found it useful to include a fourth mesh in z, the vertical distance by which the orbit itself misses the center of the target. This quantity can be used to cut events not originating from the target. The previous remarks made for p: also apply for z. Specifically, no new transformations are needed for z; also, if Y ; < 0 , the sign of z found after interpelation must be changed.
Ioo
z
ARM TWO p~55 GeV/c
ARM ONE p~17 GeV/c
80 [
w
~ ~o cr Z
40
z
© ~-O 15
,
-OIO
-005
O
005
OtO
Oi5 -oi5
-oio
-005
O
O 05
O '10
rl OI5
I00 (p- p ' ) / p
Fig. 3. Histogram of interpolation error (in %) for 1000 Monte-Carlo events, p = original lab momentum, p ' = interpolated momentum.
FAST RETRACE
OF P A R T I C L E
To conclude this section we present the linear interpolation equation, which is formulated to minimize the number of arithmetic operations. In this formulation the values of the momentum components and z at each mesh point are not used directly; instead certain linear combinations of the values at 16 neighboring mesh points are used and must be calculated in aclvance. To simplify things, let P stand for one of the three momentum components or z. After transforming an experimental event given by the four wire coordinates, we must find into which 4-dimensional rectangle the data point falls. Let Vk0 and A V k be the first mesh point and the (constant) mesh interval along the V k axis. Let I k --
v~-V~o AV~
Xk = t k - - [ t k ] ,
(13)
where [t] is the largest integer less than or equal to t. T]he four [tk] show which rectangle the data point is in while the Xk give the relative position inside the rectangle. We are assuming that the data point actually falls somewhere inside the pre-established mesh. Let P~B~6be the value of P at one of the 16 corners of the 4-dimensional rectangle; the kth subscript value is zero or one depending on whether P~a~o is on the low or high side of the rectangle along the Vk axis. The 4-dimensional linear interpolation is then given by P = C1 + C2 xl
+ ( C 3 + C 4 X l ) x 2 -~-
+ [C5 + C6 x~ +(C7 + Cs xl)x2] x3 + -'F'-{ C 9 + C l o X 1 -~- ( C l i -{- C 12 XI)X2 -~-
-+-[C13-[-C14x I + ( C I 5 + C I 6 X l ) X 2 3 X 3 } X 4 ,
(14)
where the Ci are given in table 1.
3. Speed and accuracy A typical application uses a mesh that is 5 by 5 by 2 by 2 for V~ through I/4, which amounts to only 100 mesh points and sixteen 4-dimensional rectangles. The small number of mesh points is an important consideration since much computer time is needed to calculate the momentum for each mesh point (about ¼ of a second on the CDC 7600). Also note that the space requirements for tLe mesh in the computer is small; with 16 rectangles, 1024 memory locations are needed for each spectrometer arm. The computing time needed to process data for each spectrometer arm is 0.6 ms/event (on the C D C 6600). This includes the time needed to transform the four
325
ORBITS TABLE 1
The 16 i nt e rpol a t i on c ons t a nt s used in eq. (14).
i 1 2 3 4 5 6 7 8 9 lO Il 12 13 14 15 16
G Poooo Plooo -- C1 Poloo - C1 P l l o O - (C1 + C 2 + C3) Poolo - C1 PlOlO -- (C1 + C2 + C5) P011o - (C1 + Ca + C5) 7 Pl110 - ~ C, l=l PO001 -- C1 PlOOl-(CI+C2+C9) POlOl - (C1 + C3 + Co)
Pl101-(C1+C2+C3+C4+C9+C10+Cll) POOll - (C1 + C5 + C9) P l O l l - (C1 + C2 + C5 + C6 + C9 + Clo + ClZ) Pol 11 - (C1 + C3 + C5 + C7 + C9 + C11 + C13) 15 P l l l l -- t__~lC~
wire coordinates and to interpolate all three momentum components and z. Lechanoine et al. 2) quote times of 6 and 15 ms/event (also on a CDC 6600) for two different cases that they describe. The difference in speed is not hard to find. Their use of a function of five observables consisting of sums of products of Tchebycheff polynomials is really a more complicated interpolation scheme in a five dimensional mesh of calculated orbits. (The mesh points are chosen according to certain properties of Tchebycheff polynomials.) Hence, a large number of mesh points are used in their calculation while we use only the nearest 16 points. With 100 mesh points the error introduced by the method is less than 0.5% which compares favorably with the accuracy quoted in ref. 2. A Monte-Carlo simulation was done to check the accuracy. Events chosen uniformly over a typical spectrometer setting were traced forward by numerical integration and retraced by the mesh. A histogram of the interpolation errors in both arms is shown in fig. 3. The error is over an order of magnitude smaller than the experimental resolution of the spectrometer. The biggest source of error occurs in the interpolation along the V z axis due to the rapid variation of momentum with Va for a fixed VI. This is why we generally use as many mesh peints along V 2 as along V1. Although never tried, the accuracy may be improved by rotating the lines of constant V~ so that they end up perpendicular to the central line X~ = g ( X E ) (fig. 2). Even though in
326
R.J.
our example the n u m b e r o f mesh points along V3 and V4 are as few as possible, i n t e r p o l a t i o n errors along these two axes are the smallest o f all; this is due to the already mentioned fact that the horizontal m o m e n t u m c o m p o n e n t is insensitive to V 3 and //4 while Pz and z are nearly linear in the vertical wire coordinates. 4. Conclusion To summarize: A multi-dimensional rectangular mesh is chosen in a space o f new variables whose dimensionality is equal to the n u m b e r of observables needed to define an orbit. A suitable set o f transformations between the new variables and the observables is chosen so that the mesh covers the physically meaningful orbits and linearizes at least one momentum component. Other t r a n s f o r m a t i o n s may be used to linearize other m o m e n t u m components. Once the m o m e n t u m c o m p o n e n t s are calculated at each mesh point prior to d a t a processing, the m o m e n t u m of
McKEE
each experimental event is found by t r a n s f o r m i n g the observables into the new variables and linearily interpolating in the multi-dimensional mesh. The detailed d e v e l o p m e n t of section 2 was presented to show how one might go a b o u t applying the idea to a specific apparatus. The description in section 2 might suggest that this retrace m e t h o d is strongly tied to a specific apparatus. But the basic idea - suitable t r a n s f o r m a t i o n o f observables and linear interpolation in a mesh f o r m e d by the t r a n s f o r m e d variables to find the m o m e n t u m quickly and accurately - certainly has general application.
References ~) H. L. Anderson, D. Larson, L. Myrianthopoulos, L. Dubal, C. K. Hargrove, E. P. Hincks, R. J. McKee, H. Mes, D. Kessler, and A. C. Thompson, to be published. 2) C. Lechanoine, M. Martin, and H. Wind, Nucl. Instr. and Meth. 69 (1969) 122.