A Dynamic Nonlinear Impedance Model of a Single Photon Avalanche Diode

A nonlinear LCR parallel circuit model of a single photon avalanche diode (SPAD) is derived from a Lienard-type nonlinear differential equation. The resistance and the inductance associated with avalanche multiplication (AM) are time-dependent and governed by the avalanche time constant due to the impact ionization ratio. Time dependences of current, voltage, resistance, and inductance in the model are analyzed by numerically solving the equation. During the initial generation of a Geiger-mode pulse, when the voltage reaches the breakdown voltage, the resistance diverges to limit the maximum current and the inductance reduces to give the maximum speed of the voltage variation with the avalanche time constant. In the frequency domain, avalanche impedance (AI) spectra are obtained as a ratio of voltage and current spectra calculated as Fourier transforms of the time domain signals. The AI spectra exhibit a negative-resistance and an inductance. An analytic expression for the impedance is derived and found to comprise only a quenching resistance and a stray capacitance as indicated by a constant radius of a Nyquist plot. Finally, the present model is shown to incorporate the models of impact-ionization-avalanche transit-time (IMPATT) diodes when a very small-signal limit is assumed.

A Dynamic Nonlinear Impedance Model of a Single Photon Avalanche Diode Akito Inoue , Shinzo Koyama, and Yutaka Hirose Abstract-A nonlinear LCR parallel circuit model of a single photon avalanche diode (SPAD) is derived from a Lienard-type nonlinear differential equation. The resistance and the inductance associated with avalanche multiplication (AM) are time-dependent and governed by the avalanche time constant due to the impact ionization ratio. Time dependences of current, voltage, resistance, and inductance in the model are analyzed by numerically solving the equation. During the initial generation of a Geiger-mode pulse, when the voltage reaches the breakdown voltage, the resistance diverges to limit the maximum current and the inductance reduces to give the maximum speed of the voltage variation with the avalanche time constant. In the frequency domain, avalanche impedance (AI) spectra are obtained as a ratio of voltage and current spectra calculated as Fourier transforms of the time domain signals. The AI spectra exhibit a negative-resistance and an inductance. An analytic expression for the impedance is derived and found to comprise only a quenching resistance and a stray capacitance as indicated by a constant radius of a Nyquist plot. Finally, the present model is shown to incorporate the models of impact-ionization-avalanche transit-time (IMPATT) diodes when a very small-signal limit is assumed.

I. INTRODUCTION
S INGLE-PHOTON avalanche diode (SPAD)-based CMOS image sensors (CISs) have been developed intensively in recent years as evidenced by aggressive pixel scaling [1], [2], [3], [4], by faster and more precise photon counting circuits [4], [5], [6], [7], by well-controlled quenching circuits [6], [8], [9], by incorporation of a pinned photodiode (PPD) into pixels [10], and by realization of 3-D integrated device structures [1], [4], [5], [11]. In parallel, the range of applications has been extending from time-of-flight sensors [1], [2], [3], [4], [10] and sensitive photon counting imagers [1], [2], [6], [8] to high dynamic range imagers [2], [5], [6]. Among the technological advancements, the realization of well-controlled pixel-levelquenching is most important to imager performance. end, the development of a pixel circuit model of a SPAD that represents complex carrier dynamics of a Geiger-mode pulse and, at the same time, incorporates common circuit concepts is essential. So far, mainly three types of models have been reported as follows: 1) Models consisting of discrete passive elements, switches, and transistors without intrinsic avalanching capability [12], [13], [14]. In these models, a Geigermode pulse must be input from the external supply as a voltage or current trigger pulse by controlling the switch. 2) Models in which carrier evolution is numerically calculated from basic carrier-dynamics equations [15], [16], [17]. With this method, because the number of carriers in the basic equation is treated as a continuous variable, forced carrier extinction must be introduced to simulate a successful quenching [16]. 3) Models in which carrier numbers are computed as outcomes of successive random events of impact ionizations [18], [19]. However, none of the above models describes the most pronounced physical characteristics of the avalanche-generated carriers, i.e., negative resistance and inductive effects [20], [21]. These effects would significantly affect the maximum speed and amount of avalanche carrier generation in a SPAD as extensively investigated in impact-ionization-avalanche transit time (IMPATT) diodes [21]. The roles of negative resistance and inductive effects in SPADs are largely unknown.
In this article, from a basic Lienard-type nonlinear differential equation (LNDE) [22], we develop a dynamic-(time-dependent) circuit model of a SPAD. An equivalent circuit model derived from LNDE gives a time-dependent parallel pair of a resistance and an inductance. These elements are shown to originate from the avalanche effect being inversely proportional to the avalanche time constant. From LNDE, time-dependent avalanche currents and voltages are numerically calculated with quenching resistances as parameters. It is clarified that the model rigorously specifies both the maximum voltage speed and the maximum current at a characteristic time termed the breakdown time. By taking the ratio of the voltage across and current flowing through the model SPAD in the frequency domain, avalanche impedance (AI) spectra are obtained. The AI spectra exhibit that the real and imaginary parts are, indeed, negative and inductive, respectively, in all frequency regions. In addition, AI is, then, analytically shown to be identical to a  (c) Equivalent circuit of (a) derived from (1)-(3). Conduction current (I cond ), displacement current (I disp ), and total current of SPAD (I SPAD ) flow through red, blue, and green dashed rectangular elements, respectively. C SPAD is a capacitance of MR. R A , and L A are derived as (2) and (3), respectively. C stray is a parasitic capacitance of the cathode node.
"π -phase-shifted" parallel pair of a quenching resistance and a capacitance as verified by Nyquist plots analyses. Finally, by taking a small-signal limit, the present model is shown to lead to circuit models of IMPATT diodes confirming the identity of the physical origins of the two types of devices.

II. TIME DOMAIN ANALYSES OF A CIRCUIT MODEL A. Basic Equations and an Equivalent Circuit of a SPAD
We consider a 1-D SPAD with a p-i-n structure where the i-region is a multiplication region (MR). We assume that, inside MR, an electric field due to a bias is uniform and that the width of MR is constant. A circuit model where a quenching resistor, R, is connected to the cathode of the SPAD is illustrated in Fig. 1(a). A bias voltage V 0 is above the breakdown voltage (V BD ) with a small excess voltage An anode of the SPAD is grounded. In Fig. 1(b), a band diagram of the 1-D SPAD during avalanche multiplication (AM) is illustrated. Relevant equivalent circuit elements are overlapped on the band diagram. A parallel pair of a resistor [R A (t)] and an inductor [L A (t)] are due to impact ionization inside MR as derived shortly. The capacitance C of the cathode node is a sum of the capacitances of the SPAD MR and additional stray components, i.e., C = C SPAD + C stray . We also assume that: 1) by impact ionization, carriers are generated at the center of MR and 2) carriers transport with their saturation velocities (v e for electrons and v h for holes) inside MR. The voltage across the SPAD [V (t)] is described by the following Lienard-type differential equation (in [22, eqs. (14) and (15)]): We define the coefficients ofV (t) and V (t) − V BD as, respectively, where τ Q is called as a quenching lifetime and defined in [22] as Since (4) describes the maximum avalanche generation rate when the SPAD's internal field (voltage) becomes the breakdown field (voltage), hereafter, τ Q is renamed as the avalanche time constant.
Thus, (1) represents a parallel circuit comprising R A (t), L A (t), and C = C SPAD + C stray as illustrated in Fig. 1(c) equivalent to Fig. 1(b). It should be noted that both R A (t) and L A (t) depend on τ Q indicating that these parameters directly describe the avalanche carrier dynamics. R A (t) is regarded as a resistance for current due to avalanche-generated carriers to charge/discharge C. L A (t) is considered as an inductance due to variation of the avalanche current. It is also mentioned that R A (t) and L A (t) diverge when V (t) = V BD and when V (t) = V 0 , respectively, limiting the current and its variation (speed).

B. Derivation of Current Components
Integrating (1), current conservation is derived as follows: where I C (t), I CSPAD (t), I Cstray (t), I R A (t), I L A (t), I R (t), and I 0 are current components flowing through C, C SPAD , C stray , R A (t), L A (t), R, and the current source, respectively. [I R A (t)+ I L A (t)] is the conduction current through the red dashed rectangle in Fig. 1(c) and is described as avalanche charges drifting out of MR as where n(t) is the number of electrons inside MR. The second equality of (13) is derived from [22, eq. (4)], . In addition, displacement current flowing through C SPAD is obtained as where D is an electric flux density in C SPAD . A total current flowing through the SPAD is then

C. Time Domain Analyses of Voltages and Currents
Equations (1) and (6)- (15) were solved by setting up difference equations from (1) with parameters listed in Table I. V BD and V ex are set at 27.5 V and 1.5 V, respectively, as in [22]. An initial condition is the generation of one electron-hole pair in the middle of MR. From the difference equations, V (t) was obtained at an increment of each 2.0 ps. I SPAD (t), I C (t), I R A (t), I L A (t) and I R (t) are calculated from each V (t). To reproduce a successful quenching state, n(t) is treated as a continuous variable when n(t) > 1 and is set to be zero immediately after when n(t) = 1 [16], [22]. Because the initial condition is only determined by a single charge injection, the present model is applicable to both photo-and dark current-initiated Geiger-mode pulse generation. It is also noted that, because the original model does not include carrier recombination terms, any associated effects such as after-pulse effects due to trapping and de-trapping and light emission effects due to radiative recombination, are not dealt with in the following. In Fig. 2(a), calculated V (t) with a successful (R = 75 k ) and an unsuccessful (R = 60 k ) quenching conditions are plotted as blue and red solid lines, respectively. The threshold value for successful quenching, R th , is 66 k calculated from (34) of [22]. Regardless of the quenching conditions, after the buildup period, at the breakdown time, t ≡ t BD =152 ps, V (t) drops to V BD . t BD is a period during which the avalanche-generated carriers become sizable and [V 0 − V (t)] to assume a voltage swing of V ex . t BD is typically several hundreds of picoseconds [22]. It is emphasized that, up to t ∼ 2t BD , both V (t) and I SPAD (t) are dominated by AM and insensitive to quenching conditions. We observe that the slope of V (t), i.e.,V (t), is maximum at t BD , giving the avalanche time constant τ Q , which will be observed in Fig. 3(b) shortly. Beyond V BD , V (t) drops further to reach the full voltage drop ( V Q ) nearly twice of V ex ; V Q ≃ 2V ex = 3.0V. When R = 75 k , after reaching the minimum, V (t) recovers to V 0 with a time constant of RC, much longer than τ Q . With R = 60 k , quenching fails since n(t) does not fall below unity according to the unsuccessful quenching condition. As a result, V (t) oscillates around V BD with decreasing envelope and approaches a steady state condition.
In Fig. 2(b) and (c), the time dependence of current flowing through each circuit element of Fig. 1(c) is plotted. I cond (t) (green curves) is seen to be supplied by I C (t) (orange curves) as indicated by their opposite signs. They reach the maximum and the minimum values, respectively, at t BD = 152 ps. Thus, discharge of the capacitances C(= C SPAD +C Stray ) is the source of I cond (t) consistent with [13] and [15]. For recharging C, I C (t) (orange curves) are supplied by I R (t) (blue curves) as indicated by their opposite signs. It is noted that the results in Fig. 2(b), a successful quenching, are only possible when I cond (t) is forced to be zero before V (t) recovers above V BD . Signals of I SPAD (t) are plotted by red lines which comprise pulse-shaped conduction current, displacement current, and long-tailed recharge current. Reduction in I SPAD (t) compared with I cond (t) is due to I disp (t).With R = 75 k , after I cond (t) drops to 0, I R A (t) (yellow dashed line) and I L A (t) (light blue dashed line) are due to induced currents by charges flowing out of C canceling each other. The effect is regarded as reciprocal of the Ramo-Shockley theorem [19].
With R = 60 k , failed quenching signals of repeated AM pulses with gradually decreasing height are well reproduced. Both I cond (t) and I SPAD (t) approach the steady state with a constant value, V ex /R. In Fig. 3, R A (t) and L A (t) calculated from (2) and (3) are plotted. In both graphs, orange and green lines correspond to, respectively, the successful (R = 75 k ) and the unsuccessful (R = 60 k ) quenching conditions. It is noted that R A (t) is negative when V (t) > V BD . R A (t) diverges when V (t) = V BD consistent with the fact that V BD is an equilibrium point of the system giving dn(t)/dt = 0 and therefore, d I cond (t)/dt = 0. Thus, the divergence of R A (t) sets the maximum I cond (t) and therefore the maximum avalanche charges. It should also be noted that the time of the R A (t) divergence specifies the effective avalanche time, t BD as indicated by a vertical dotted  Fig. 3(a). On the other hand, L A (t) decreases from infinite at t = 0 to τ Q R at t = t BD and to the minimum of τ Q R/2 when V (t) reaches the minimum point and then, gradually increases during the recharge phase. We confirm that L A (t BD )/R = τ Q , the maximum speed of V (t) as observed in Fig. 2(a). The above equality is analytically verified by [22, eq. (36)] and by (22) in Section IV-B derived as the time constant of IMPATT diodes which operate under the identical bias condition. It is emphasized that the analyses with prespecified t BD and R th as above are impossible in previous SPAD models where artificial avalanche charges are externally supplied [12], [13], [14].

D. Time Domain Analyses of Energy
By integrating both sides of (5) multiplied by V (t), the energy contained in each circuit element is calculated. The results are shown in Fig. 4(a) (R = 75 k ) and Fig. 4(b) (R = 60 k ), where E C (t) and E R (t), respectively, denote energies contained in C, and R. E cond (t) and E I (t), are energies due to I cond (t) and I 0 , respectively. E C (t), [E R (t) − E I (t)], and E cond (t) are plotted as red, blue, and green solid lines, respectively. When quenching is successful [see Fig. 4(a)], as with I C (t) and I cond (t), at t = t BD , E cond (t) increases stepwise in exchange for a sudden decrease of E C (t). Then, E C (t) recovers to the initial value supplied by [E R (t)− E I (t)], whereas E cond (t) stays constant being in the idling state. The total energy change of the SPAD, i.e., E loss = |E R (∞) − E I (∞)| = |E cond (∞)|, is the same amount as the maximum decrease of E C (t), calculated as follows: Equation (16) is twice as much as that described in [13, Section III-H]. The difference comes from the estimation of voltage drops, i.e., V Q = V ex in [13] and V Q = 2V ex in this work. Equation (16) is interpreted as necessary work (energy) for a total amount of generated charge (2C V ex ) transferred against the breakdown voltage. When quenching is failed [see Fig. 4(b)], E cond (t) increases stepwise as AMs occur repeatedly. E C (t) never recovers to the initial value.

III. FREQUENCY ANALYSES OF AI
We calculate AI spectra by taking the ratio of the Fourier transforms of V (t) and I SPAD (t) as  Using V FFT ( f ), I FFT ( f ), and (17), AI spectra based on FFT, Z FFT ( f ), are calculated. The results are shown as Nyquist plots together with the corresponding real and imaginary parts in Fig. 6. It is noted that the real and imaginary parts of Z FFT ( f ) are always negative and positive (inductive), respectively, consistent with the characteristics of AM and the inertial effect of the current [20], [21]. Correspondingly, the Nyquist plots are located entirely in the second quadrant of the complex plain. In addition, the plots form semicircles of constant radius of R/2 centered at (−R/2, 0).

IV. DISCUSSIONS A. Explicit Form of AI
The Nyquist plots with constant radii (see Fig. 6) are verified as follows. From (13)-(15) Taking Fourier transforms of both sides of (18) leads to (20) which is independent of device parameters. In Fig. 6, the Nyquist plots calculated from (20) are plotted as solid lines for R = 75 k and 60 k . It is clearly seen that agreement between Z ( f ) and Z FFT ( f ) is nearly perfect. This is due to the current conservation (5) where the current flowing through R and C Stray is in the opposite phase to I SPAD (t). It is also noted that the negative sign of (20) indicates a phase shift by π from the phase of an ordinary parallel RC pair. Finally, a Nyquist plot with R th = 66 k (black dashed line) specifies the critical boundary outside of which quenching is successful. Such clear specification of R th with a given storage capacitance, C, and a quenching time, τ Q [see (4)], was unavailable with previous models where quenching resistances were selected based essentially on a rule of thumb, e.g., a few tens of k or larger up to 1 M with a large scattering in values [12], [13], [14], [15].
V. CONCLUSION A nonlinear equivalent circuit model of the avalanche breakdown in a SPAD is derived [see Fig. 1(c)] from an LNDE (1). Explicit expressions for a negative resistance R A (t) (2) and an inductance L A (t) (3) are derived. Both R A (t) and L A (t) are shown to be governed by the avalanche time constant τ Q . Time-dependent currents and voltage in each element are calculated and confirmed to be conserved (see Fig. 2). The maximum current and speed of the Geiger-mode pulse are found to be limited by the divergence of R A (t) (see Figs. 2 and 3) at the breakdown time, t = t BD . The frequency spectra I FFT ( f ) and V FFT ( f ) calculated from the time-domain calculations (see Fig. 5) are found to critically depend on the quenching conditions. The Nyquist plots of the AI, Z FFT ( f ), (see Fig. 6) clearly exhibit the negative resistance and an inductance. An analytical form of AI, i.e., Z ( f ), is derived as (20) and confirmed to be identical to Z FFT ( f ). Z ( f ) is also found to be expressed only by the external components of the model SPAD, i.e., R and C stray (20). Finally, the present SPAD model is shown to incorporate those of IMPATT diodes when a small-signal limit is taken [(21) and (22)]. The present method enables one to design SPAD pixel circuits with a priori knowledge of avalanche time constants, maximum of generated charges, and quenching conditions unavailable with other models.