A track reconstruction method for silicon microstrip detectors

NUCLEAR INSTRUMENTS 8 METHODS IN PHVSICS RESEARCH ELSEWIER

Sectm

Nuclear Instruments

and Methods in Physics Research

A

A 399 (1997) 27-34

A track reconstruction method for silicon microstrip detectors Tim McMahon * Langston

Universit.v. P.O. Box 217. Langston. Received

OK 73050. USA

22 April 1997

Abstract A method for reconstructing the position of tracks traversing a silicon microstrip detector is presented. The method attempts a maximal use of the information in the charge collected at individual strips. The position resolution is determined for various strip pitch, signal to noise and incident track angle values. Application to pixel detectors is also discussed.

1. Introduction Silicon microstrip and pixel detectors have continued to be a preferred choice for experiments requiring high precision tracking due to their superior spatial resolving capabilities [l-7]. We outline here an algorithm for reconstructing track position which attempts to make a maximal use of the information in the charge collected at the strips in order to obtain the highest possible resolution. Various methods for reconstructing track position have been presented to date [8,9]. The method presented here is a maximum likelihood approach which finds the track position that has the highest probability of generating the distribution of charge seen at the strips. This approach has a number of advantages. First, the same algorithm is used for all track incident angles. Second, overlapping hits from two tracks can be reconstructed and improve the two hit resolution. Third, clusters of strips can be fit simultaneously for track position and incident angle and the angle used to aid in global event reconstruction. * E-mail: [email protected]. 0168-9002/97/$17.00 0 1997 Elsevier Science B.V. All rights reserved PIISO168-9002(97)00921-2

Fourth, given sufficient statistics the diffusion and energy loss parameters could be extracted from the data and be used to improve the reconstruction. Finally, the method is easily generalized from the one-dimensional microstrip case to a two-dimensional pixel. The paper is arranged as follows. Section 2 discusses the Monte-Carlo simulation used to generate the data for the study. In that section we also compare the results of our model for charge transport with CERN test beam data. Section 3 discusses the algorithm developed here for track reconstruction on single track and double track events. Section 4 gives the results of using the algorithm for the track resolution obtained for 25 and 50 pm strip pitch detectors.

2. Silicon strip simulation For the energy deposition of an ionizing particle in the silicon strip detector we use the model from Bichsel [lo] and start with the single collision differential cross section a(E) and take 45 GeVlc pions as the simulated particle. The cross section for multiple

T. McMahonlNucl.

28

Instr. and Meth. in Phys. Res. A 399 (1997) 27-34

0.0150 L

strips

(4 NT -

o.oooo

4

0

50

100 energy

200

150

loss

250

300

/( Fig. 2. Diagram of simulated detector. The origin is placed at the track exit. Strip regions traversed by the track are labeled T strips and those not traversed labeled NT.

A (eV)

0))

where a(A)*’ = a(A) for an energy loss A. For a thickness t of silicon the energy loss probability is given by a Poisson weighting of the n collision convolution functions

0

400

200

energy

600

loss

800

(c)

0.0006 -

0

1000

2000

energy Fig. 1. Energy loss function 5 pm thick slice of silicon.

fa

loss for (a)

J

3000

A (eV)

1pm, (b) 2vm and (c)

collisions is given by the n-fold convolution single collision cross section[ IO],

of the

A

u(a)*”

=

0

@)o*(“-‘)(A

f~(t, A) = 2 y n=O

,(A)*‘, .

A (eV)

0.0008 -

A

1000

- E) dE,

(1)

where m, the average number of collisions in an absorber of thickness t, is given by m = MO where MO, the number of collisions per unit length, is taken as 38400cm-’ for 45 GeVlc pions[lO]. To model a 300 pm thick detector we calculated the convolutions using Eq. (1) up to n = 1500 where each of the ~(a)*~ functions were discretized in 1500 steps in In(E). Figs. 1 show f A for thicknesses of 1, 2 and 5 pm. We have neglected the small loss in the energy the incident particle undergoes after traversing the thin slices t when constructing the CT(A)*~ functions. A positive feature of the convolution method is that the average energy loss per unit thickness will be independent of the simulation step length chosen. This property is useful in our case since the simulation is modeling the energy deposited in strips from thicknesses of 0.5 to 300 pm. In the simulation, tracking of secondary delta rays and multiple scattering is not done. A diagram of the detector and coordinate system used is shown in Fig. 2. We define the strip region as the area underneath the strip where the width is

T McMahonj

Nucl

Instr. and Meth.

in Phvs. Rex A 399

f 1997) 27-34

indicated with the dotted line. Charge drifts vertically upward under the influence of a bias voltage. Strip regions which are physically traversed by the particle are labeled as T strips as indicated and non-traversed regions are labeled as NT strips. From diffusion and capacitive cross-talk, signals may appear on NT strips from the track. Simulated tracks are stepped through the detector with steps sizes of at maximum 5 urn. Energy loss during the step is calculated once and is spread uniformly along the length of the track step segment. The distribution of the charge collected from this track segment at the strips is given by an analytic formula derived as follows. The charge density at the surface from a track segment length dl with incident angle 0 at a depth .v is given by a Gaussian distribution

0 -100

-80

de(x)= dx

1 2 sin 0

discontinuous

equation

&

(A

=

b-

fi -(x

A=m

+ vtan 0).

(51

The & signs in Eq. (4) are the solutions for x > 0 and x < 0, respectively[l2]. The equation has a singularity at 6 = 0’ and for this angle an equation can be found by applying I’Hospital’s rule. Fig. 3 shows the charge profile dQ/dr for a 300um thick slice of silicon at various track angles. For large angles the diffusion has an asymmetric smearing of charge to the left of the cluster. The total charge collected at a strip is found by integrating Eq. (4) and yields the

0

20

40

for various track angles.

(eC7”an’t.y’D( ‘f I

s])

-erf[A]) . 1)

(6)

where the i signs in Eq. (6) are the solutions x > 0 and x < 0, respectively, and A given by Eq. (5 ). The discontinuity at .X= 0 must be removed when the range includes this point. In that case where the endpoints are x_ = x < 0 and x+ = .r > 0 the total charge is Ptotal = (P(-e)

with

-20 (pm)

Fig. 3. Surface charge density dQ!dr

-erf coefficient. The drift vefrom taking the electric as constant. Integrating for a track starting at 3: the charge density at the

-40 x

Q(x)

where D is the diffusion locity I’ is made constant field in the bulk detector Eq. (3) over all dl pieces and ending at y = 0 yields surface

-60

- ecx-

11+ (!&+)

- Qt +c)l

(7)

with t: a small amount away from 0. The equation also has a singularity for incident tracks with (3 = 0” and should use an angle approaching 0 = 0”. These equations are defined with the end point of the track segment at the surface or y = 0. For track segments not ending at J’ = 0 the equations are evaluated twice, once at each end point -vendand _vstart,and the difference between the two taken. If a magnetic field is present then the charge drifts at a Lorentz angle 0~ away from vertical. A solution in this case is found by applying to Eqs. (4) and (6) the transformations y + _V/‘cosHL, 0-0’=tan-‘(cos0~(tanl~r+tant9)).

(8)

30

T McMahon/Nucl.

Instr. and Meth. in Phys. Res. A 399 (1997) 27-34

The dQ/d.x and Q values obtained after the y and f3 translation are finally transformed by

(cos;;;;s”‘) de/&, cosOLcos0’ Q- ( case)Q.

dQ1d.x +

(9)

Eqs. (4) and (6) are easily generalized to a rectangular surface area, a pixel, under the assumption that diffusion is independent in the x and z direction. The equations are evaluated separately for the two dimensions. The final two-dimensional probability function is given by multiplying the two together for a joint probability. We have fit CERN silicon microstrip test beam data [I 1J with Eq. (4) as shown by the solid lines in Fig. 4. The dashed lines are from a numerically integrated solution of Eq. (3) where the electric field is not taken as a constant but varies linearly in the silicon from (vb - Vd)/t to (vb + &)/t with Vb the bias voltage and vd the depletion voltage. For the fit the drift mobility of holes was set to ~1 = 450 cm*/Vs and only the diffusion coefficient D, a horizontal offset and an overall normalization were left as free parameters. When usmg Eq. (4) the drift velocity was fixed to 3.21 x 1O6cm/s for the zero B-field data and 1.93 x lo6 cm/s for the B = 1.68 T data. Fig. 4(a) is data with no magnetic field and the fit gives D = 15.9 & 0.8 cm2/s using the average electric field formula, Eq. (4), and D = 19.4 f 1.1 cm2/s using the linear electric field formula. For the 1.68 T data the Lorentz angle was added as a free parameter. Fig. 4(b) is data with a B = 1.68 T magnetic field and yields D = 15.7 k 1.9 cm*/s and Lorentz angle 8L = 2.50 III 0.07’. The linear E-field formula gives D = 19.0 f 2.5 cm2/s and Lorentz angle iJL = 2.86 i 0.08”. Events were simulated for 25 and 50 urn strip pitch detectors and track incident angles from 0” to 65”. D, the diffusion coefficient, was set to 35 cm*/s and 11, the drift velocity, set to 1.35 x 1O6cm/s. The signal to noise ratio (S/N) is defined as the total average charge deposited by the track divided by the noise charge in a single strip. Gaussian noise is added to each strip with 3500eV for S/N = 25, 1750 eV for SIN = 50 and OeV for S/N = 00.

x f/d

-15

-10

-5

0

5

10

15

x (vd

Fig. 4. Fit of surface charge density function to CERN test beam data with constant (solid) and linear (dashed) electric field charge transport functions for (a) no magnetic field and (b) 1.68 T magnetic field.

3. Track reconstruction The reconstruction of track position for the MonteCarlo events from above proceeds by placing tracks at horizontal positions in the range of the strips which have collected charge. For each track the energy deposited along the track is unfolded from the strip charge taking into account charge diffusion. Probabilities for the energies deposited to result from the thickness of the region traversed are used to calculate a likelihood for that particular horizontal position.

31

T. McMahon I Nucl. Instr. and Muth. in Phys. Rrs. A 399 (1997) 27-34

All simulated tracks were required to be reconstructed. The strip cluster finding starts with the strip that has the largest signal and is in proximity to the generated track. Strips to the left and right are added to the cluster until the strip falls below threshold. On both sides, the first strip that was below threshold was also finally added to the cluster. The threshold cut was set to 2800eV for all S/N values used. The selected strips are put into an nctuster element vector Q. The center of gravity (COG) of the clusters were calculated and serve as a point of reference of the resolution of the Monte-Carlo simulation and gauge the improvement from the charge profile (QP) fitting algorithm. 3.1. Singltbtruck

nr matrix Vi/ that will now map E back to an (nT + 3 ) element vector Q’ by Q’; =

C

V,j ET.

j=l

This study did not include capacitive cross-talk. It can, when determined for a particular detector, be included as well to the T and V matrices. The likelihood for the values in the vector E is computed from the noise convoluted fA functions for the two left most and two right most strips in E as

reconstruction

(12)

The QP reconstruction method finds the horizontal position of a track that gives the maximum likelihood of generating the signals Q observed in the cluster. In this study the incident angle of the fitted track is fixed to the generated value and is a reasonable constraint since typically the track angle is given from global tracking that includes hit positions from additional silicon microstrips and/or other detector devices. A 100 pm range of track positions around the cluster is scanned. A track placed at a particular trial position will physically traverse nT number of the strip regions. For this track an nr x nr matrix T is calculated which maps how deposited charge along a track segment in a strip region is transported and distributed amongst the nr surface strips including the effects of diffusion. Eqs. (6) and (7) are used to determine this mapping. The rows of T,j run over surface strip number and the columns over the traversed strip region. An II~ element vector E of the predicted deposited charge along the track is then calculated from

(10) where only the nr elements of Q are used and joffset depends on the offset of the exit of the track relative to the end of the cluster. The track lengths of the vector E are also determined and put to the vector L. Next, the T matrix is expanded by adding the two surface strips to the left and one to the right for an (nr + 3) x

Inner strips are not included to the likelihood since they give marginal information about the horizontal position of the track and add noise. For the nontraversed (NT) strips the likelihood is calculated from the two left and one right strips in Q’ with a Gaussian probability

with gnoise the single-strip noise. The number of T and NT strips the likelihoods are calculated for can vary depending upon the placement of the candidate track relative to the cluster. Therefore, an average In C over the strips used is calculated and then used as the final criteria for selecting the best candidate track from the scan. Some track clusters present special cases, for example, frcluster may be less than nr, but, they are easily addressed with some additional logic. 3.2. Overlapping

trucks

Reconstructing clusters which are from two proximate or overlapping tracks proceeds in a similar way as the single-track case including some extensions. Two tracks are horizontally placed at trial locations and a likelihood calculated. Again a matrix T is calculated using Eqs. (6) and (7) that maps how a unit charge deposited by both tracks is transported and distributed over the surface strips. For regions traversed by both tracks the unit charge is partitioned between

32

T McMahon/

Nucl. Instr. and Meth. in Phys. Rex A 399 (1997) 27-34

the two track segments in proportion to the length of the tracks in the region. As is done in Eq. (10) a vector E is calculated from multiplying Ton Q where the nri&r2 elements used from Q are the union of the traversed strips of tracks 1 and 2. Vectors El and E2 for the two tracks are extracted from the vector E from Lli E1i = (Lli + LZi) Eh

25

(a) 0

(14)

where L2i in Eq. (14) will be zero if track 2 is not traversing region i. The vector E2 is made analogously as El in Eq. (14). The Q’ vector is then calculated as ml

Qli = C

x4XxX

m2

Vl,

El, + C

j=l

V&j E2j,

o-““““““““““““.“”

0

(15)

where Vl and V2 are the T matrices extended to include extra surrounding strips. The In C values are then calculated analogously as Eqs. (12) and (13) before where terms from both tracks will be present. For this study only single-track events were generated. The simulation and track reconstruction was run on 50 processors concurrently in the PVM environment on an IBM-SP2 at the Maui High Performance Computing Center. Parallelism was achieved by parsing events among separate processors and was scalable to an arbitrary number of processors.

4.1. Track position In Fig. 5 we show the position resolutions obtained from the COG and QP methods for a signal to noise level of 25 and 50 for a 25 urn pitch detector. The resolutions for the QP method at the angle closest to 0” were 1.3 urn and 0.7 pm for the 25 and 50 S/N values, respectively, and can be compared with the corresponding numbers of the COG method of 2.2 pm and 1.4 pm. The values were obtained by fitting the track residuals histogram with a Gaussian for a given angle. All simulated tracks are required to be reconstructed and entered into the residual histogram. No goodness of fit criteria were used in cutting tracks from the sample and poorly reconstructed tracks will contribute to outliers in the Gaussian fit. The percentage of these tracks for the QP algorithm has been estimated

20

30

track

40

angle

* 0

3

e

2 2

t

50

(degrees)

10 g d

X

0

e*

5-

0 0 0

&p

4. Results

10

incident

j=l

0

&ZX

XX

XX

xx

x

QXX

xc$ ,‘,,,~‘~~,~‘~,~,‘,~~,‘,~~~‘~~ 0

x

XXX

10

incident

20

track

30

angle

40

50

I

(degrees)

Fig. 5. Resolution of 25 pm pitch detector with signal to noise of (a) SiN = 25 and (b) S/N = 50 for the COG (diamond) and QP (cross) methods.

by counting tracks that are beyond 3 sigma from the mean. For the S/N = 25 case, 2-5% of the tracks are beyond 3 sigma for angles O-30”. Above 30” the percentage increases to 37% at 50”. For S/N = 50 the equivalent numbers are l-2% for O-30” and increase to 8% at 50”. For the simulation with no noise l-2% of tracks are beyond 3 sigma for all angles. Fig. 6 shows the resolution obtained if no strip noise is added to the simulated tracks. Error bars are statistical only. Fig. 7 shows the results for a 50pm pitch detector. The results obtained here are from our

T. McMahon

I Nucl. Instr. and Med

33

in Phys. Res. A 399 (1997) 27-34

(4 t

20

t

--z 3

15 i

“.”

0

10 incident

20 track

30 angle

Fig. 6. Resolution of 25 km pitch detector (diamond) and QP (cross).

40

0-----L-0

50

20 incident

(degrees) with no noise. COG

20

60

40 track

angle

track

angle

(degrees)

(b) t

software implementation of the algorithm have room for improvement.

and may

4.2. Truck skewing

As can be seen in Fig. 3 tracks with large angles will have their charge spread asymmetrically from the effects of diffusion. One would anticipate that the COG method would systematically shift the reconstructed position to the left because of this. We have studied this effect and have found that this initial view is complicated by the threshold cut used to find clusters. This threshold cut can in fact result in systematic shifts towards the right since the diffused charge that falls below the threshold, is cut and the remaining cluster has a net movement of the COG towards the right. The resulting track position skewing has an angular dependence and is also dependent on the threshold cut applied. In the analysis we performed the threshold cut was set to 2800eV. For the COG method the track skewing for the 25 pm pitch case was less than 0.1 urn at 0” and increased linearly to approximately f0.5 urn at 40’. Track skewing was also seen in the QP method with a similar angular dependence, but, at a level of half the magnitude of the COG numbers.

20

0

incident

40

60 (degrees)

Fig. 7. Resolution of 50 pm pitch detector with signal to noise of (a) S/N = 25 and (b) S/N = 50 for the COG (diamond) and QP (cross) methods.

Acknowledgements I thank Prof. George Kalbfleisch at the Univ. of Oklahoma for useful discussions and pointing out the Bichsel work and Eric Smith who made a check of the test beam data. The undergraduate assistants at Langston University John Timberlake, Renolda Grant, Luhua Lin, Billy Gaston, and Eddie Walker worked on running simulations and rendering plots. Thanks to Wayne Trail of Kinema Research and Software for reading this. Dr. Frank Gilfeather at the Univ. of New Mexico and Prof. Dennis Judd at Prairie View A&M

34

T. McMahon 1 Nucl. Instr. and Meth. in Phys. Rex A 399 (1997) 27-34

Univ. helped to provide computer time at the Maui High Performance Computer Center. This work was supported in part by a grant from the Dept. of Energy, grant # DE-FG02-95ER40923. Support also came through the use of the Maui High Performance Computing Center, sponsored by the Phillips Laboratory, Air Force Materiel Command, USAF, under cooperative agreement number F2960193-2-0001. The views and conclusions contained in this document are those of the author and should not be interpreted as necessarily representing the official policies or endorsements, either expressed or implied, of Phillips Laboratory, the U.S. Government, The University of New Mexico, or the Maui High Performance Computing Center.

References [I] M. Acciarri et al., Nucl. Instr. and Meth. A 351 (1994) 300. [2] R. Lipton et al., DO Silicon Tracker Technical Design Report, 5117194.

[3] D. Amidei et al., Nucl. Instr. and Meth. A 350 (1994) 73. [4] J. Blocki et al., NucI. Instr. and Meth. A 342 (1994) 269. [5] W.W. Armstrong et al., Atlas Technical Proposal, CERNILHCC194-43, 12/l 5194. [6] G.L. Bayatian et al., CMS Technical Proposal, CERN/LHCC 94-38, 12/15/94. [7] M. Campbell et al., Nucl. Instr. and Meth. A 342 (1994) 52. [S] R. Turchetta, Nucl. Instr. and Meth. A 335 (1993) 44. [9] N. Colonna and E. Lisi, Nucl. Instr. and Meth. A 334 (1993) 551. [lo] Hans Bichsel, Rev. Modem Phys. 60 (3) (1988). [l I] E. Belau et al., Nucl. Instr. and Meth. 214 (1983) 253. [12] To derive the equations use was made of Mathematics, Wolfram Research Inc. and the Handbook of Mathematical Functions, Milton Abramowitz and Irene Stegun, National Bureau of Standards.