Theoretical and simulation studies of characteristics of a Compton light source

Compton scattering of a laser beam with a relativistic electron beam has been used to generate intense, highly polarized and nearly monoenergetic x-ray or gamma-ray beams at many facilities. The ability to predict the spatial, spectral and temporal characteristics of a Compton gamma-ray beam is crucial for the optimization of the operation of a Compton light source as well as for the applications utilizing the Compton beam. In this paper, we present two approaches, one based upon analytical calculations and the other based upon Monte Carlo simulations, to study the Compton scattering process for various electron and laser beam parameters as well as different gamma-beam collimation conditions. These approaches have been successfully applied to characterize Compton gamma-ray beams, after being benchmarked against experimental results at the High Intensity Gamma-ray Source (HIGS) facility at Duke University.


I. INTRODUCTION
Compton scattering of a laser beam with a relativistic electron beam has been successfully used to generate intense, highly polarized and nearly monoenergetic x-ray or gamma-ray beams with a tunable energy at many facilities [1][2][3]. These unique Compton photon beams have been used in a wide range of basic and application research fields from nuclear physics to astrophysics, from medical research to homeland security and industrial applications [1].
The ability to predict the spectral, spatial and temporal characteristics of a Compton gamma-ray beam is crucial for the optimization of the gamma-ray beam production as well as for research applications utilizing the beam. While the theory of particle-particle (or electronphoton) Compton scattering, which is equivalent to the scattering between a monoenergetic electron beam and a monoenergetic laser beam with zero transverse sizes, is well documented in literature [4][5][6], there remains a need to fully understand the characteristics of the gamma-ray beam produced by Compton scattering of a laser beam and an electron beam with specific spatial and energy distributions, i.e., the beam-beam scattering.
Study of beam-beam Compton scattering has been recently reported in [7,8]. However, the algorithms used in these works are based upon the Thomson scattering cross section, i.e., an elastic scattering of electromagnetic radiation by a charged particle without the recoil effect. For scattering of a high energy electron beam and a laser beam, the recoil of the electron must be taken into account. The Compton scattering cross section has been used to study characteristics of Compton gammaray beams by Duke scientists in 1990's [9,10]. However, the effects of incoming beam parameters, and the effects of gamma-beam collimation were not fully taken into account.
In this paper, we present two different methods, a semi-analytical calculation and a Monte Carlo simulation, to study the Compton scattering process of a polarized (or unpolarized) laser beam with an unpolarized electron beam in the linear Compton scattering regime. Using these two methods, we are able to characterize a Compton gamma-ray beam with various laser and electron beam parameters, arbitrary collision angles, and different gamma-beam collimation conditions. This paper is organized as follows. In Section II, we first review the calculation of the Compton scattered photon energy for an arbitrary collision angle, and then introduce the scattering cross section in a Lorentz invariant form. Based upon this cross section, the spatial and spectral distributions as well as the polarization of a Compton gamma-ray beam are investigated in particleparticle scattering cases. In Section III, we discuss the beam-beam Compton scattering by considering effects of the incoming beam parameters as well as the effect of the gamma-ray beam collimation. Two methods, a semi-analytical calculation and a Monte Carlo simulation, are then presented. Based upon the algorithms of these methods, two computing codes, a numerical integration code and a Monte Carlo simulation code, have been developed at Duke University. The benchmarking results and applications of these two codes are presented in Section IV. The summary is given in Section V.

II. PARTICLE-PARTICLE SCATTERING
A. Scattered photon energy A review of the calculation of scattered photon energies in the particle-particle scattering case is in order. and a photon in a lab frame coordinate system (xe, ye, ze) in which the electron is incident along the ze-direction. The incident photon is propagating along the direction given by the polar angle θi and azimuthal angle φi. The collision occurs at the origin of the coordinate system. After the scattering, the scattered photon propagates in the direction given by the polar angle θ f and azimuthal angle φ f . θp is the angle between the momenta of incident and scattered photons, k and k ′ . The electron after scattering is not shown in the figure.
electron and a photon in a laboratory frame coordinate system (x e , y e , z e ) in which the incident electron with a momentum p is moving along the z e -direction. The incident photon with a momentum k ( is the Planck constant) is propagated along the direction with angles (θ i , φ i ). The collision occurs at the origin of the coordinate system. After the collision, the photon with a momentum k ′ is scattered into the direction of (θ f , φ f ).
According to the conservation of the 4-momenta before and after scattering, we can have where p = (E e /c, p) and k = (E p /c, k) are the 4momenta of the electron and photon before the scattering, respectively; p ′ = (E ′ e /c, p ′ ) and k ′ = (E g /c, k ′ ) are their 4-momenta after the scattering; E e and E p are the energies of the electron and photon before the scattering; E ′ e and E g are their energies after the scattering; and c is the speed of light. Squaring both sides of Eq. (1) and following some simple manipulations, we can obtain the scattered photon energy as follows, where β = v/c is the speed of the incident electron relative to the speed of light, and θ p is the angle between the momenta of the incident and scattered photons (Fig. 1). For a head-on collision, θ i = π and θ p = π − θ f , Eq. (2) can be simplified to Clearly, given the energies of the incident electron and photon, E e and E p , the scattered photon energy E g only depends on the scattering angle θ f , independent of the azimuth angle φ f . The relation between the scattered photon energy E g and scattering angle θ f is demonstrated in Fig. 2. In this figure, the scattered photon energies E g are indicated by the quantities associated with the concentric circles in the observation plane, and the scattering angles θ f are represented by the radii R of the circles, i.e, θ f = R/L, where L = 60 meters is the distance between the collision point and the observation plane. We can see that the scattered photons with higher energies are concentrated around the center (θ f = 0), while lower energy photons are distributed away from the center. Such a relation, in principle, allows the formation of a scattered photon beam with a small energy-spread using a simple geometrical collimation technique.
For a small scattering angle (θ f ≪ 1) and an ultrarelativistic electron (γ ≫ 1), Eq. (3) can be simplified to where γ = E e /(mc 2 ) is the Lorentz factor of the electron and mc 2 is its rest energy. When the photon is scattered into the backward direction of the incident photon (i.e., θ f = 0, sometimes called backscattering), the scattered photon energy will reach the maximum value given by Neglecting the recoil effect, i.e., 4γ 2 E p /E e ≪ 1, Eq. (5) can be reduced to the result given by the relativistic Thomson scattering theory [8] We can see that the incident photon energy E p is boosted by a factor of approximately 4γ 2 after the backscattering. Therefore, the Compton scattering of photons with relativistic electrons can be used to produce high energy photons, i.e., gamma-ray photons. Under a set of conditions θ i ≈ π and θ f ≈ 0, the uncertainties of the scattered photon energy E g due to the uncertainties of the variables (E e , E p , θ f and θ i ) in Eq. (2) can be estimated [10,11]. For example, the relative uncertainty of the scattered photon energy ∆E g /E g due to the uncertainty of the electron beam energy ∆E e /E e is given by taking the derivative of Eq. (2) with respect to E e , i.e., Contributions to ∆E g /E g associated with other variables are summarized in Table I.

Lorentz invariant form
The general problem concerning the collision is to find the probabilities of final states for a given initial state of the system, i.e., the scattering cross section. Using Quantum Electrodynamics (QED) theory, the Compton scattering cross section in the Lorentz invariant form has been calculated in [4,12,13], and the result for unpolarized electrons scattering with polarized photons is given by where r e is the classical electron radius; φ f is the azimuthal angle of the scattered photon; ξ 1,2,3 and ξ ′

1,2,3
are Stokes parameters describing the incident and scattered photon polarizations in their respective coordinate systems; X and Y are the Lorentz invariant variables defined as follows where s and u are the Mandelstam variables [4] given by X and Y satisfy the inequalities [4] Since the scattering cross section of Eq. (8) is expressed in the Lorentz invariants, it can easily be expressed in terms of the collision parameters defined in any specific frame of reference.

Polarization description in lab frame
In the laboratory frame, three right-hand coordinate systems are used in Eq. (8) to describe the motion and polarization of the incident electron (x e , y e , z e ), the incident photon (x,ỹ,z), and the scattered photon (x ′ ,ỹ ′ ,z ′ ) (Fig. 3). The coordinate system (x e , y e , z e ) is fixed in the lab frame, and its z e -axis is along the incident direction of the electron. (x,ỹ,z) and (x ′ ,ỹ ′ ,z ′ ) are the local coor- (xe, ye, ze) is the coordinate system for the incident electron ( p) moving along the ze-axis direction. For the head-on collision, the incident photon ( k) comes along the negative ze-axis, and the scattered photon ( k ′ ) is moving along the direction given by the polar angle θ f and azimuthal angle φ f . The momentum vectors k and k ′ form the scattering plane. (x,ỹ,z) is a right-hand coordinate system attached to the scattering plane. Thezaxis is along the direction of k;x-axis is perpendicular to the scatter plane; andỹ-axis is in the scattering plane. (x ′ ,ỹ ′ ,z ′ ) is another right-hand coordinate system attached to the scattering plane. Thez ′ -axis is along the direction of k ′ ;x ′ -axis is the same as thex-axis; andỹ ′ -axis lies in the scattering plane.
dinate systems attached to the scattering plane formed by the momenta of the incident and scattered photons, k and k ′ . For (x,ỹ,z), thex-axis is perpendicular to the scattering plane; theỹ-andz-axes are in the scattering plane with thez-axis along the direction of k. For (x ′ ,ỹ ′ ,z ′ ), thex ′ -axis is the same as thex-axis for the incident photon, perpendicular to the scattering plane; and thez ′ -axis is along the direction of k ′ .
The Stokes parameters ξ The polarization of the photon is always defined in its local coordinate system with its momentum being one of the axes. For Compton scattering described by Eq. (8), these local coordinate systems (x,ỹ,z) and (x ′ ,ỹ ′ ,z ′ ) are different for different scattering planes. However, for the cases that the photons and electrons collide nearly head on to produce high-energy photons with small scattering angles, it becomes possible to conveniently express in an approximate manner the polarization of the incident and scattered photons using a fixed coordinate system, for example, the lab-frame electron coordinate system (x e , y e , z e ).
Let us consider the incident photon with itsz-axis approximately parallel to the negative z e -axis. The Stokes parameter of the incident photon can be related to the degrees of polarization defined in the fixed electron coordinate system through the following equations [5,14], where P t and P c are the degree of linear and circular polarizations of the incident photon defined in the coordinate system (x e , y e , z e ), respectively; τ is the azimuthal angle of the linear polarization P t with respect to the x e -axis; and φ f is the azimuthal angle of the scattering plane. For Compton scattering involving an ultra-relativistic electron, scattered photons are concentrated in a small scattering angle (θ f < 1/γ). For these high-energy photons with small scattering angles, theirz ′ -axes are approximately parallel to the z e -axis. Neglecting the polar angle (i.e. θ f ≪ 1), the Stokes parameters of the scattered photon can be expressed approximately using a set of Stokes parameters defined in the fixed electron coordinate system as [14], whereξ ′ 1,2,3 are the Stokes parameters defined in the coordinate system (x e , y e , z e ).

C. Spatial and energy distributions of scattered photons
Based upon Eqs. (8), (12) and (13), we can calculate the spatial and energy distributions of a gamma-ray beam produced by Compton scattering of a monoenergetic electron and laser beams with zero transverse beam sizes, i.e., the particle-particle scattering.
Let us consider Compton scattering of an unpolarized electron and a polarized laser photon without regard to their polarizations after the scattering. The differential cross section is obtained by setting ξ ′ 1,2,3 to zero in Eq. (8) and multiplying the result by a factor of two for the summation over the polarizations of the scattered photons [4]. Thus, the differential cross section is given by [11] dσ The total cross section can be obtained by integrating Eq. (14) with respect to Y and φ f , Note that the Stokes parameter ξ 3 depends on φ f ; however, after integration over φ f the dependence vanishes.
Neglecting the recoil effect (X ≪ 1), we can have which is just the classical Thomson cross section.

Spatial distribution
For a head-on collision (θ i = π) in a laboratory frame, according to Eq. (9) the Lorentz invariant quantities X and Y are given by and Substituting dY into Eq. (14), the angular differential cross section is given by where dΩ = sin θ f dθ f dφ f and ξ 3 has been expressed in terms of P t (Eq. (12)). From Eq. (19), we can see that the differential cross section depends on the azimuthal angle φ f of the scattered photon through the term P t cos(2τ − 2φ f ). For a circularly polarized or unpolarized incident photon beam (P t = 0), this dependency vanishes. Therefore, the distribution of scattered photons is azimuthally symmetric. However, for a linearly polarized incident photon beam (P t = 0), the differential cross section is azimuthally modulated, and the gamma photon distribution is azimuthally asymmetric. Figs. 4 and 5 illustrate the spatial distributions of Compton gamma photons at a location 60 meters downstream from the collision point for both circularly and linearly polarized incident photon beams. In these figures we can also see that the distribution of scattered photons peaks sharply along the direction of the incident electron beam. This demonstrates that the gamma-ray photons produced by Compton scattering of a relativistic electron beam and a laser beam are mostly scattered into the electron beam direction within a narrow cone.

Energy distribution
For a head-on collision in the laboratory frame, it can be shown that Thus, Substituting dY in Eq. (14) and integrating the result with respect to the azimuth angle φ f , we can obtain the energy distribution of scattered photons as follows The energy spectrum calculated using Eq. (22) is shown in Fig. 6. The spectrum has a high energy cutoff edge which is determined by the incident electron and photon energies according to Eq. (5). In Fig. 6, we can see the spectral intensity has a maximum value at the scattering angle θ f = 0, and a minimum value around the scattering angle θ f = 1/γ. The ratio between them is about 2 when the recoil effect is negligible. This will be shown in the next section.
Note that the energy spectrum shown in Fig. 6 is for a Compton gamma-ray beam without collimation. However, if the gamma-ray beam is collimated by a round aperture with a radius of R and distance L from the collision point, the energy spectrum will have a low energy cutoff edge, and its value can be calculated using Eq. (4) with θ f = R/L.

Observations for a small recoil effect
For a small recoil effect (X ≪ 1), we can approximate Eqs. (19) and (22) to draw several useful conclusions. For convenience, we first define Using the inequality Eq. (11), it can be found that approximately (with a negligible recoil effect, X ≪ 1), Thus, the maximum and minimum spectral flux of the Compton gamma-ray beam are given by and The ratio between them is which is shown in Fig. 6. When θ f = 0, we can have Substituting Y in Eq. (23), we have f (Y ) ≈ 1/2. Thus, the spectral flux has a maximum value around the scattering angle θ f = 0. When θ f = 1/γ, we can have Substituting Y into Eq. (23), we have f (Y ) ≈ 1/4. Therefore, the spectral flux has a minimum value around the scattering angle θ f = 1/γ. These results are illustrated in Fig. 6. Expressed in terms of the total scattering cross section of Eq. (16), the fraction of scattered photons in the energy range [E max This is a simple formula which can be used to estimate the portion of the total gamma-ray flux with a desirable energy spread ∆E max g after collimation.
For a circularly polarized or unpolarized incident photon beam, according to Eq. (19), it can also be calculated that the angular intensity of scattered gamma-ray photons at the scattering angle θ f = 1/γ is about 1/8 of the maximum intensity at the scattering angle θ f = 0, i.e., In addition, integrating Eq. (14) over the entire solid angle of the cone with a half-opening angle of 1/γ, i.e., integrating Y over the range of X(1−X/2) Y X(1− X) and φ f over the range from 0 to 2π, we can have Comparing Eq. (33) to the total cross section of Eq. (16), we can conclude that about half of the total gamma-ray photons are scattered into the 1/γ cone. This can be explained by considering the Compton scattering in the electron rest frame. In this frame, the Compton scattering process is just like "dipole" radiation: the gammaray photons are scattered in all the directions, a half of the gamma photons is scattered into the forward direction, and the other half into the backward direction. When transformed to the laboratory frame, the gammaray photon scattered into the forward direction in the rest frame will be concentrated in the 1/γ cone in the laboratory frame.

D. Polarization of scattered photons
For polarized photons scattering with unpolarized electrons without regard to the final electron polarization, the cross section is given by Eq. (8). Substituting ξ 1,2,3 and ξ ′ 1,2,3 using Eqs. (12) and (13), and assuming the linear polarization of the incident photon beam is along the x e -axis, i.e., τ = 0, we can get It should be noted that the Stokes parametersξ ′ 1,2,3 describe the polarization of the scattered photon selected by a detector, not the polarization of the photon itself [4].
In order to distinguish them from the detected Stokes parametersξ ′ 1,2,3 , we denote the Stokes parameters of the scattered photon itself by ξ f 1,2,3 . According to the rules presented in section 65 of [4], ξ f 1,2,3 are given by Integrating Eq. (34) over the azimuthal angle φ f gives where Therefore, the averaged Stokes parameters of the scattered photons over the angle φ f are given by ξ f i = Φ i / Φ 0 , which depend on the incident photon polarization and variables X and Y .
For example, for 100% horizontally polarized (P t = 1, P c = 0, τ = 0) incident photons scattering with unpolarized electrons, the average Stokes parameters of the Clearly, the scattered photons retain the polarization of the incident photons. ξ f 3 as a function of the scattered photon energy is shown in Fig. 7 for 800 nm laser photons head-on colliding with 500 MeV electrons. It can be seen that the average Stokes parameter ξ f 3 of scattered gamma-ray photons is almost equal to 1 around the maximum scattered photon energy as in this case the recoil effect is negligible. It means the scattered gamma-ray photons with the maximum energy are almost 100% horizontally polarized.

III. BEAM-BEAM SCATTERING
In the previous section we discussed the spatial and spectral distributions of a gamma-ray beam produced by Compton scattering of monoenergetic electron and laser beams with zero transverse beam sizes, i.e., particleparticle scattering. However, in the reality, the incoming electron and laser beams have finite spatial and energy distributions, which will change the distributions of the scattered gamma-ray beam. Therefore, there remains a need to understand the characteristics of a Compton gamma-ray beam produced by scattering of a laser beam and an electron beam with specific spatial and energy distributions, i.e., the beam-beam scattering.
In this section, we discuss the beam-beam Compton scattering process. First, we derive a simple formula to Two coordinate systems are defined to describe electron and laser beams: the first coordinate system (x, y, z) is the electron-beam coordinate system in which the electron beam is moving along the z-axis direction; the (x l , y l , z l ) system is the laser-beam coordinate system in which the laser beam propagates in the negative z l -axis direction. The coordinate systems (x, y, z) and (x l , y l , z l ) share the same origin.
calculate the total flux of the Compton gamma-ray beam. Then, we present two methods, a semi-analytical calculation and a Monte Carlo simulation, to study the spatial and spectral distributions of the gamma-ray beam. Based upon these methods, two computing codes, a numerical integration code and a Monte Carlo simulation code have been developed. These two codes have been benchmarked against the experimental results at High Intensity Gamma-ray Source (HIγS) facility at Duke University.
A. Geometry of beam-beam scattering Figure 8 shows Compton scattering of a pulsed electron beam and a pulsed laser beam in a laboratory frame. Two coordinate systems are used: (x, y, z) for the electronbeam moving along the z-direction; the (x l , y l , z l ) for the laser-beam propagating in the negative z l -direction. These two coordinate systems share a common origin. The time t = 0 is chosen for the instant when the centers of the electron beam and laser pulse arrive at the origin. The definition of these two coordinate systems allows the study of the Compton scattering process with an arbitrary collision angle, i.e, the angle between z-axis and negative z l -axis. For a head-on collision, the collision angle equals π. In this case, the electron and laser coordinate systems coincide.
In these coordinate systems, the electron and laser beams with Gaussian distributions in their phase spaces can be described by their respective intensity functions as follows [9] f e (x, y, z, x ′ , y ′ , p, t)= 1 (2π) 3 p is the momentum of an electron, and p 0 is the centroid momentum of the electron beam; x ′ and y ′ are the angular divergences of the electron beam in the x-and ydirections, respectively; α x,y , β x,y and γ x,y are Twiss parameters of the electron beam; σ p , σ z and ε x,y are the electron beam momentum spread, RMS bunch length, and transverse emittance, respectively; k and λ are the wavenumber and wavelength of a laser photon, and k 0 is the centroid wavenumber of the laser beam; β 0 , σ k and σ l are the Rayleigh range, the RMS energy spread and bunch length of the laser beam. Note that the waist of the laser beam is assumed to be at the origin of both coordinate systems.

B. Total flux
The number of collisions occurring during a time dt and inside a phase space volume d 3 p d 3 k dV is given by [5] dN ( r, p, k, t) = σ tot ( p, k)c(1 − β · k/| k|)n e ( r, p, t) where σ tot ( p, k) is the total Compton scattering cross section which is determined by the momenta of the incident electron and laser photon, p and k; β = v e /c is the relative velocity of the incident electron; n e ( r, p, t) = N e f e ( r, p, t) and n p ( r, k, t) = N p f p ( r, k, t), where f e ( r, p, t) and f p ( r, k, t) are the phase space intensity functions of electron beam and laser pulse, and N e and N p are the total numbers of electrons and laser photons in their respective pulses.
To calculate the total number of scattered gamma-ray photons produced by collision, Eq. (41) needs to be integrated for the entire phase space and the collision time, i.e., where θ i is the collision angle between the incident electron and laser photon. Assuming collisions occur at the waists of both beams (α x = α y = 0, σ w = λβ 0 /(4π)), the spatial and momentum phase space in the density functions can be separated, i.e, f e ( r, p, t) = f e ( r, t)f e ( p) and f p ( r, k, t) = f p ( r, t)f p ( k). Since the cross section σ tot ( p, k) only depends on p and k, we can have where is the single-collision luminosity defined as the number of scattering events produced per unit scattering cross section, which has dimensions of 1/area [15]. For a headon collision (θ i = π) of a relativistic electron (β ≈ 1) and a photon, the single-collision luminosity can be simplified to Thus, Eq. (43) can be rewritten in a simple form where σ tot is the total Compton scattering cross section averaged over the momenta of the incident electrons and photons. Neglecting the energy spread of the electrons and photons, σ tot can be approximated by σ tot of Eq. (15), which can be further simplified to the classical Thomson cross section if the recoil effect is negligible. If the beam-beam collision rate is f 0 , the gamma-ray flux is given by C. Spatial and energy distributions: semi-analytical calculation To obtain the spatial and energy distributions of a Compton gamma-ray beam, the differential cross section should be used instead of the total cross section in Eq. (42). In addition, two constraints need to be imposed during the integration of Eq. (42) [9,10]. First, let us consider the geometric constraint, which assures the gamma-ray photon generated at the location r can reach the location r d shown in Fig. 9. In terms of the position vector, this constraint is given by where k ′ represents the momentum of the gamma-ray photon; r = (x, y, z) denotes the location of the collision; and r d = (x d , y d , z d ) denotes the location where the scattered gamma-ray photon is detected. Due to the finite spatial distribution and angular divergence of the electron beam, a gamma-ray photon reaching the location r d can be scattered from an electron at different collision points with different angular divergences. The constraint of Eq. (48) projected in the x-z and y-z planes is given by Here, θ x and θ y are the projections of the scattering angle θ f in the x-z and y-z planes, i.e., θ x = θ f cos φ f , θ y = θ f sin φ f and θ 2 f = θ 2 x + θ 2 y , where θ f and φ f are the angles defined in the electron coordinate system (x e , y e , z e ) in which the electron is incident along the z edirection (Fig. 9). x ′ and y ′ are the angular divergences of the incident electron, i.e., the angles between the electron momentum and z-axis. L is the distance between the collision point and the detection plane (or the collimation plane). Note that a far field detection (or collimation) has been assumed, i.e., L ≫ | r| and L ≈ | r d |.
The second constraint is the energy conservation. Due to the finite energy spread of the electron beam, the gamma-ray photon with an energy of E g can be produced by electrons with various energies and scattering angles. Mathematically, this constraint is given by Imposing the geometric and energy constraints in Eq. (42), the spatial and energy distributions of a Compton gamma-ray beam can be obtained by integrating all the individual scattering events, i.e., where dΩ d = dx d dy d /L 2 , and dσ/dΩ is the differential Compton scattering cross section. Note that a head-on collision between electron and laser beams has been assumed, and the density function f e ( r, p, t) has been replaced with f e (x, y, z, x ′ , y ′ , p, t) of Eq. (40) under the approximation p z ≈ p for a relativistic electron beam. In addition, the integration · · · f p ( r, k, t)d 3 k is replaced with · · · f p (x, y, z, k, t)dk, where f p (x, y, z, k, t) is defined in Eq. (40). Integrations over dk x and dk y have been carried out since the differential cross section has a very weak dependency on k x and k y for a relativistic electron beam.
Assuming head-on collisions for each individual scattering event (θ i = π and dσ/dΩ is given by Eq. (19)), neglecting the angular divergences of the laser beam and replacing x ′ and y ′ with θ x and θ y , we can integrate Eq. (52) over dV, dt and dp to yield the following result (see Appendix A), where and σ Ee is the RMS energy spread of the electron beam. In a storage ring, the vertical emittance of the electron beam is typically much smaller than the horizontal emittance. For a Compton scattering occurring at a location with similar horizontal and vertical beta functions (β x ∼ β y ), the vertical divergence of the electron beam can be neglected. In addition, the photon energy spread of a laser beam is small, and its impact can also be neglected in many practical cases. Under these circumstances, the cross section term in Eq. (53) has a weak dependence on θ y (≈ y d /L) and k (≈ k 0 ). With the assumption of an unpolarized or circularly polarized laser beam, Eq. (53) can be simplified further after integrating θ y and k: where The integrations with respect to k, θ y and θ x in Eq. (53) or θ x in Eq. (55) must be carried out numerically. For this purpose, a numerical integration Compton scattering code (CCSC) in the C++ computing language has been developed to evaluate the integrals of Eqs. (53) and (55). With the detailed spatial and energy distributions of the Compton gamma-ray beam dN (E g , x d , y d )/(dE g dx d dy d ), the energy spectrum of the gamma-ray beam collimated by a round aperture with a radius of R can be easily obtained by integrating dN (E g , x d , y d )/(dE g dx d dy d ) over the variables x d and y d for the entire opening aperture, i.e., x 2 d + y 2 d R 2 . The transverse misalignment effect of the collimator on the gamma-ray beam distributions can be introduced by replacing x d and y d with x d + ∆x and y d + ∆y in Eq. (53) or Eq. (55), where ∆x and ∆y are the collimator offset errors in the horizontal and vertical directions, respectively.

D. Spatial and energy distributions: Monte Carlo simulation
In the previous section, we have derived an analytical formula to study the spatial and energy distributions of a Compton gamma-ray beam. However, to simplify the calculation several approximations have been made: head-on collisions for each individual scattering event, a negligible angular divergence of the laser beam, and far field collimation.
A completely different approach to study Compton scattering process is to use a Monte Carlo simulation. With this numerical technique, effects that cannot be easily included in an analytical method can be properly accounted for. For example, using a Monte Carlo simulation we can study the scattering process for an arbitrary collision angle. With this motivation, we developed a Monte Carlo Compton scattering code. In following, the algorithm of this code is presented.

Simulation setup
At the beginning of the collision, both the electron and laser pulses are located some distance away from the origin (Fig. 8), and two pulse centers arrive at the origin at the same time (t = 0). The collision duration is divided into a number of time steps, and the time step number represents the time in the simulation.
Due to a large number of electrons in the bunch, it is not practical to track each electron in the simulation. Therefore, the electron bunch is divided into a number of macro particles (for example, 10 6 ) which are tracked in the simulation.
The phase space coordinates of each macro particle are sampled at time t = 0. For an electron beam with Gaussian distributions in phase space, the coordinates are sampled according to the electron beam Twiss parameters as follows [16,17] x(0) = 2u 1 ε x β x cos φ 1 , y(0) = 2u 2 ε y β y cos φ 2 , y ′ (0) = − 2u 2 ε y /β y (α y cos φ 2 + sin φ 2 ), where u 1,2 are random numbers generated using an exponential distribution with a unit mean parameter (i.e., e −u1,2 ), r 1,2 are random numbers generated according to a Gaussian distribution with a zero mean and unit standard deviation, and φ 1,2 are uniformly distributed random numbers between 0 and 2π. The coordinates of macro particles at any other time (t = 0) can then be obtained by transforming the coordinates given by Eq. (56). The Compton scattering is simulated according to the local intensity and momentum of the laser beam at the collision point. The intensity of the laser beam at the collision point (x, y, z) in the electron-beam coordinate system can be calculated according to Eq. (40) using the laser-beam coordinates (x l , y l , z l ) transformed from (x, y, z). The momentum directionk of the photon at the collision point (x, y, z) can be calculated from the point of view of electromagnetic wave of the photon beam. For a Gaussian laser beam, its propagation phase ψ(x l , y l , z l ) in the laser-beam coordinate system is given by [16,18] ψ(x l , y l , z l ) = −ik l z l − ik l z l x 2 l + y 2 l 2(β 2 0 + z 2 l ) ; (57) the wavevector (the momentum of photon k l ) is given by k l = ▽ψ(x l , y l , z l ). Thus, where The unit vectork l expressed in the electron-beam coordinate system gives the momentum direction of the laser photon in this coordinate system.

Simulation procedures
At each time step, the Compton scattering process is simulated for each macro particle. The simulation proceeds in two stages. In the first stage, the scattering probability is calculated using the local intensity and momentum of the laser beam. According to this probability, the scattering event is sampled. If the scattering happens, a gamma-ray photon will be generated, and the simulation proceeds to the next stage. In the second stage, the energy and scattering angles (including the polar and azimuthal angles) of the gamma-ray photon are sampled according to the differential Compton scattering cross section. The detailed simulation procedures for these two stages are presented as follows.

First stage: scattering event
Since the energy and scattering angles of the gammaray photon are not the concern at this stage, the total scattering cross section is used to calculate the scattering probability. According to Eq. (41), the scattering probability P ( r, p, k, t) in the time step ∆t for the macro particle at the collision point (x, y, z) is given by where n p (x, y, z, k, t) and k are the local density and wavevector of the photon beam, respectively; σ tot ( p, k) is the total scattering cross section given by Eq. (15).
According to the probability P ( r, p, k, t), the scattering event is sampled using the rejection method as follows [19,20]: first, a random number r 3 is uniformly generated in the range from 0 to 1; if r 3 ≤ P ( r, p, k, t), Compton scattering happens; otherwise the scattering does not happen, and the above sampling process is repeated for the next macro particle.

Second stage: scattered photon energy and direction
When a Compton scattering event happens, a gammaray photon is generated. The simulation proceeds to the next stage to determine the energy and scattering angles of the gamma-ray photon. For convenience, the sampling probability for generating gamma-ray photon parameters is calculated in the electron-rest frame coordinate system (x ′ e , y ′ e , z ′ e ) in which the electron is at rest and the laser photon is propagated along the z ′ e -axis direction. Since the momenta of macro particles and laser photons have been expressed in the electron-beam coordinate system (x, y, z) in the lab frame, we need to transform the momenta to those defined in the electron-rest frame coordinate system (x ′ e , y ′ e , z ′ e ). After transformations, the sampling probability for generating the scattered gamma-ray photon energy and direction will be calculated as follows.
In the electron-rest frame coordinate system (x ′ e , y ′ e , z ′ e ), according to Eq. (2) the scattered photon energy is given by where θ ′ is the scattering angle between the momenta of the scattered and incident photons; E ′ g and E ′ p are the energies of the scattered and incident photons, and E ′ g is in the range of In the electron-rest frame coordinate system, we can simplify the Lorentz invariant quantities X and Y of Eq. (14) to X = 2E ′ p /mc 2 and Y = 2E ′ g /mc 2 . As a result, the differential cross section is given by where τ ′ is the azimuthal angle of the linear polarization direction of the incident photon beam defined in the system (x ′ e , y ′ e , z ′ e ), and φ ′ is the azimuthal angle of the scattered photon. Note that the quantity P t , the degree of linear polarization of the incident photon beam, is invariant under Lorentz transformations.
The scattered photon energy E ′ g and the azimuthal angle φ ′ are sampled according to the differential cross section Eq. (63). Since Eq. (63) depends on both E ′ g and φ ′ , the composition and rejection sampling method [19,20] is used to sample these two variables. To sample the scattered gamma-ray photon energy E ′ g , Eq. (63) needs to be summed over the azimuthal angle φ ′ and written as where and 0 ≤ f (E ′ g ) ≤ 1 for any E ′ g . Now, the scattered gamma-ray photon energy E ′ g can be sampled according to f (E ′ g ) as follows: first, a uniform random number E ′ g is generated in the range given by Eq. (62), and r 4 in the range from 0 to 1; if r 4 ≤ f (E ′ g ), E ′ g is accepted, otherwise the above sampling process is repeated until E ′ g is accepted. If E ′ g is accepted, the scattering angle θ ′ can be calculated using Eq. (61).
After the scattered gamma-ray photon energy E ′ g is determined, the azimuthal φ ′ angle is sampled according to The stairs plot represents the spectrum simulated using the code MCCMPT, the dash line represents the spectrum calculated using the code CCSC, and the circles represent the one using the code CAIN2.35. The electron beam energy and RMS energy spread are 400 MeV and 0.2%, respectively. The electron beam horizontal emittance is 10 nm-rad, and the vertical emittance is neglected. The laser wavelength is 600 nm with negligible photon beam energy spread. The gamma-ray beam is collimated by an aperture with a radius of 12 mm located 60 meters downstream from the collision point.
After obtaining the gamma-ray photon energy E ′ g , and the angles θ ′ and φ ′ in the electron-rest frame coordinate system, we need to transform these parameters to those in the lab-frame coordinate system. In the meantime, the momentum of the scattered electron is also computed. This electron can still interact with the laser photon in following time steps, which allows to correctly model the multiple scattering process between the electrons and laser photons.

IV. BENCHMARK AND APPLICATIONS OF COMPTON SCATTERING CODES
Based upon the algorithms discussed in Section III, we have developed two computer codes using the C++ programming language: the numerical integration Compton scattering code CCSC and the Monte Carlo Compton scattering code MCCMPT. Below, we briefly discuss the benchmark and applications of these two codes.

A. Energy Distribution
Our Compton scattering computer codes MCCMPT 11: Comparison between the measured and calculated energy spectra of a Compton gamma-ray beam. The solid line represents the calculated spectrum using the CCSC code, and the circles represent the measured gamma-beam energy distribution after removing the escape peaks and Compton plateau using a spectrum unfolding technique. The gammaray beam is produced by Compton scattering of a 466 MeV electron beam and a 790 nm laser beam at the HIγS facility. The RMS energy spread of the electron beam is 0.1%, and horizontal and vertical emittance are 7.8 and 1.0 nm-rad, respectively. The collimator with an aperture radius of 12.7 mm is placed 60 meters downstream from the collision point.
for International Linear Collider [16]. The energy spectra of Compton gamma-ray beams generated using these three codes are shown in Fig. 10. We can see that these three codes can produce very close results. In terms of computing time, the codes CCSC, MCCMPT and CAIN2.35 took about 10 minutes, 150 minutes and 1200 minutes to generate these spectra using a singlecore machine, respectively. Compared to the multipurpose beam-beam colliding code CAIN2.35, the dedicated Compton scattering codes CCSC and MCCMPT are much faster and easy to use.
At the HIγS facility, the Compton gamma-ray beam is usually measured using a high-purity Germanium (HPGe) detector. Due to the non-ideal response of the detector, the measured spectrum has a structure of a full energy peak, a single and double escape peaks, and a Compton plateau. To unfold the measured energy spectrum, a novel end-to-end spectrum reconstruction method has been recently developed [21]. The comparison of the measured gamma spectrum and calculated spectrum using the CCSC code is shown in Fig. 11. A very good agreement between them is observed.
Using the Monte Carlo simulation code, we can study the Compton scattering process with an arbitrary collision angle. The simulated spectra using MCCMPT are compared to those using CAIN2.35 in Fig. 12. Again, very good agreements are observed. It is clearly shown that the gamma-ray beam produced by a head-on collision of an electron and a laser beams has the highest energy and flux. With a 90 degree collision angle, the maximum energy of the gamma-ray beam is only half of that for a head-on collision.
The energy spread of a Compton gamma-ray beam is mainly determined by the degree of the collimation of the gamma beam, energy spread and angular divergence of the electron beam [21]. The contributions of these parameters to the gamma-ray beam energy spread are summarized in Table I. In some literature [22,23], a simple quadratic sum of individual contributions was used to estimate the energy spread of the Compton scattering gamma-ray beam. Since the electron beam angular divergence and the gamma-beam collimation introduce non-Gaussian broadening effects on the gamma-beam spectrum [21], causing the spectrum to have a long energy tail (Figs. 10 and 11), the energy spread of the gamma-ray beam cannot be given simply by the quadrature sum of different broadening mechanisms. The realistic gammaray beam energy spread needs to be calculated from its energy spectrum, which can be done using either the numerical integration code CCSC, or a Monte Carlo simulation code, MCCMPT or CAIN2.35. Figure 13 shows the spatial distribution of a gammaray beam simulated by the MCCMPT code for circularly and linear polarized incoming laser beams. For comparison, the measured spatial distributions of gamma-ray beams using the recently developed gamma-ray imaging system at HIγS facility [24] are also shown in Fig. 13. It can be seen that for a circularly polarized incoming laser beam, the distribution is azimuthally symmetric; for a linearly polarized incoming laser beam, the gamma-ray beam distribution is asymmetric, and is "pinched" along the direction of the laser beam polarization.

B. Spatial distribution
More applications of using CCSC and MCCMPT codes to study characteristics of Compton gamma-ray beams can be found in [21,25,26].

V. SUMMARY
To study characteristics of a gamma-ray beam produced by Compton scattering of an electron beam and a laser beam, we have developed two algorithms: one based upon an analytical calculation and the other using a Monte Carlo simulation. According to these algorithms, two computer codes, a numerical integration code (CCSC) and a Monte Carlo simulation code (MCCMPT), have been developed at Duke University. These codes have been extensively benchmarked against a beam-beam colliding code CAIN2.35 developed at KEK and measurement results at the High Intensity Gamma-ray Source (HIγS) facility at Duke University. Using these two codes, we are able to characterize Compton gamma-ray beams with various electron and laser beam parameters, arbitrary collision angles, and different gamma-beam collimation conditions.
In this work, the nonlinear Compton scattering process is not considered, and the polarization of the electron beam is not taken into account. Although the polarization of the gamma-ray beam has been calculated in Section II, this calculation is limited to the particleparticle scattering case. Further studies will be carried out to address these issues. where ξ x = 1 + (α x − β x L ) 2 + 2kβ x ε x β 0 , ζ x = 1 + 2kβ x ε x β 0 , σ θx = ε x ξ x β x ζ x , ξ y = 1 + (α y − β y L ) 2 + 2kβ y ε y β 0 , ζ y = 1 + 2kβ y ε y β 0 , σ θy = ε y ξ y β y ζ y , Next, we need to integrate the electron beam momentum dp. It is convenient to change the momentum p to the scaled electron beam energy variableγ = E e /(mc 2 ), and rewrite the delta-function δ(Ē g − E g ) as where is the root of with the condition of 0 ≤ θ f ≤ 4γ 2 E p E g (1 + γ 2 θ 2 f ) where