Secondary photons are emitted during the avalanche process of a single photon avalanche detector (SPAD) and they contribute to the internal and external cross-talk. In this example, we demonstrate how to calculate the transmission function between the secondary light source location inside SPAD and the measuring microscope objective. This transmission is a correction factor that, when combined with the measurements of the far-field secondary emission spectrum, allows us to calculate the secondary photon production spectrum at the source inside SPAD. We then demonstrate how to scale the far-field simulated power to include the correct source photon production spectrum. In the end, we demonstrate how to perform an external light absorption simulation in the same SPAD.

## Overview

Understand the simulation workflow and key results

### Step 1: Simulate the correction factor using STACK to find initial parameters

In this example, we are considering the following air/SiO
_{
2
}
/Si structure, where the light is generated in the Si layer.

The glass SiO
_{
2
}
layer of a SPAD’s structure heavily influences the correction factor. The reflections between the Si \(\rightarrow\) SiO
_{
2
}
and SiO
_{
2
}
\(\rightarrow\) air are wavelength- and angle-dependent. The light from the source is attenuated in the Si bulk before reaching the SiO
_{
2
}
layer. A thick SiO
_{
2
}
layer introduces an interference pattern that is dependent on the thickness. Since the true thickness of SiO
_{
2
}
is not known in this case, this step focuses on finding the thickness by matching the measured interference pattern [1] to the simulated one.

Please note, we use this step to fit the value of glass thickness and reproduce the measured interference pattern since we do not have a reliable thickness value for the proprietary device used in this example. You can skip this step if you already know all device parameters.

### Step 2: Simulate the correction factor using 3D FDTD to obtain the photon production curve

The thickness found in step 1 is used to perform a simulation of the 3D FDTD SPAD structure. This results in the transmission function, photon production curve, and emission vs angle and frequency. The true photon production curve is found by taking the ratio of the measured spectrum and the transmission function. This photon production spectrum is then used to scale the simulated far-field power and obtain the secondary emission as a function of angle, which matches the measurement when integrated over angles that cover the microscope objective.

### Step 3: Calculate absorption vs. external source angle

Using the parameters created in Step 1, we can simulate design implications for the choice of SiO
_{
2
}
thickness. This step calculates the short circuit current due to external light absorption as a function of plane wave source angle, excluding the avalanche multiplication gain.

## Run and Results

Instructions for running the model and discussion of key results

### Step 1: Simulate the correction factor using STACK to find initial parameters

- Open [[stack_purcell_sweep.lsf]].
- Manually set the location of the highest two measured peaks desired to be matched (Lines 12 and 13).
- Run [[stack_purcell_sweep.lsf]].

Light is produced from an electron avalanche in a high field region (the star point in the above image). The light travels through the silicon bulk and SiO
_{
2
}
layer before being emitted into the atmosphere. Light can be lost through reflection between layers or absorbed in the silicon bulk. The effect of a finite numerical aperture of the microscope objective needs to be considered in the far field. The STACK simulation also considers the interference effects from a thick layer of SiO
_{
2
}
. However, the actual thickness of the SiO
_{
2
}
layer is not well known but can be taken as an input parameter of the simulation.

Below we can see the measurements of the emitted light that were performed in Ref [1]; this can be used as a target result for fitting the oxide thickness:

This is the uncorrected measured emission spectrum from an FBK SiPM device [1]. The transmission through the optical stack shown in the previous figure is not yet corrected for.

The two most prominent peaks in the published data are located at 701 nm and 859 nm. The STACK sweep simulation produces the transmission for various SiO
_{
2
}
thicknesses. Finding peaks in the simulated data is accomplished using the
findpeaks
script command on the simulation result for each SiO
_{
2
}
thickness:

This figure shows the location of the highest two peaks as a function of SiO
_{
2
}
thickness. We can estimate the SiO
_{
2
}
thickness of 1.34 um because it matches best with the published data.

We can confirm this by inspecting the relative error between the measured and simulated two highest peaks:

This plot demonstrates SiO
_{
2
}
thickness vs relative error of the highest peak and the second highest peak of simulation to measurement. This method shows that the spectral data matches well at a range of 1.34 um to 1.36 um.

### Step 2: Simulate the correction factor using 3D FDTD to obtain the photon production curve

- Open fdtd_dipole_3D.fsp.
- Open fdtd_dipole_3D.lsf and set the correct SiO2 thickness discovered in Step 1.
- Run fdtd_dipole_3D.lsf.

The STACK simulation provides an initial assessment of the SiO2 thickness, but the actual thickness used for the 3D FDTD simulation will likely need to be further tweaked due to small differences in results between STACK and FDTD. Running STACK first, if the geometry allows it, still provides a lot of simulation time saving in narrowing down the parameter space. In this example, the final SiO2 thickness was set to 1.33 um to match the interference pattern.

The transmission function, the source photon production spectrum, and the emission angles are shown in figures below respectively:

Transmission curve from the 3D FDTD simulation. Another interpretation is that the y-axis is the probability of a photon leaving the avalanche region and reaching the microscope objective.

Photon production spectrum over 550nm to 1000nm. The photon production spectrum is calculated as a ratio of the measured spectrum and the FDTD transmission spectrum. The oscillations in the uncorrected measured secondary emission spectrum, shown in a figure in Step 1, have been removed. The photon production spectrum also shows larger values than the measured spectrum due to Si absorption. The photon production spectrum is then used to scale the far-field simulated power:

Integrating this angle-resolved power over the polar angle theta, and assuming symmetry in azimuthal angle phi, will result in the measured secondary emission spectrum.

### Step 3: Calculate Absorption Vs. External Source Angle

- Open fdtd_external_wave.fsp.
- Open fdtd_external_wave.lsf.
- Set the SiO2 thickness, and wavelength desired. Then run the script.

Information about the light entering the SPAD device can now be gathered. An important parameter for characterizing SiPM devices is the short circuit current for light entering at angles measured from the normal to the surface. The results of such a simulation with 2D FDTD are shown here:

Short circuit current vs light incident at an angle for 700 nm. The angle is measured normal to the surface of the SiO
_{
2
}
and the result is averaged over different polarizations. The multiplication gain is not considered in the FDTD-simulated absorption current.

## Important Model Settings

Description of important objects and settings used in this model

### Selecting Photon Source Depth

There are several important parameters to be considered when constructing a simulation model of a SPAD. The depth in which the photons are emitted are dependent on the pn junction location inside the device, which controls the location of the avalanche region. For the FBK device utilized here, the depth of the largest electric field is approximately at 0.145 um below the Si/SiO2 interface. This gives the location of the source dipole for emitted light studies. The derivation of these parameters can be found in Ref. [2].

### Photon Source Coherence

To properly simulate incoherent light produced at the source, the isotropic light source with random polarization must be simulated with multiple simulations each using a dipole oriented along a different axis. For the 3D FDTD simulation, this requires three simulations: one in the X-dir, Y-dir, and Z-dir. The results from these simulations are averaged to obtain the final incoherent solution. The external light simulations are also performed separately, as the light from the source is not coherent. Both S and P polarized light are considered in separate simulations and averaged to produce the final angle vs short circuit current result. However, in practice, the light source used in an actual measurement may not be entirely unpolarized. Care should be taken when matching the results produced through this simulation to actual measurements.

### Choice of Simulation Domain

To ensure simulations to be completed within a reasonable amount of time, the 3D FDTD simulation domain must be kept reasonably small. In this example, one method is to calculate the angle of total internal reflection and use that to define the required size of the monitor above the SiO
_{
2
}
. Or the simulation domain surface area can be verified by inspecting the real portion of the Poynting vector in that monitor. As shown below, using an example of a Z-oriented dipole simulation, all the key features are collected and the Poynting vector magnitude drops an order of magnitude near the boundaries of the simulation domain.

Magnitude of the real part of the Poynting vector through the monitor surface for the z-oriented 3D FDTD simulation.

## Updating the Model With Your Parameters

Instructions for updating the model based on your device parameters

### Optical Layer Stack

In this example the optical layer stack on top of the active Si layer consists of a single layer of SiO2 with approximately known thickness. To find the correct thickness we perform a sweep with the fast 1D STACK solver, since the simplicity of geometry allows it. If the users have a different optical layer stack they need to modify device models in fsp files accordingly, as well as make any required modifications to the fdtd post-processing lsf scripts. If all dimensions in the optical stack are known there is no need to sweep the unknown parameters to find the best fit to the measured interference pattern and step 1 can be skipped. On the other hand, if there are unknown parameters in the optical stack and the optical stack is patterned it may not be possible to do a fast sweep with the 1D STACK solver, but a full FDTD simulation sweep may need to be used.

### Si Active Layer

pn junction in Si is formed by doping and this junction determines the location of the depletion region and high electric fields required for avalanche. For the dipole location in step 1 and 2 we assume the location of the highest field, while for the absorption box thickness in step 3 we assume the depletion region width. These two parameters will depend on the doping profile and the electrical model of the device.

## Taking the Model Further

Information and tips for users that want to further customize the model

### Internal Crosstalk Simulation

Step 2 shows how to simulate the angle-resolved far field power due to secondary emission from a single SPAD. This is responsible for the external crosstalk (e.g. if there are two SPAD array chips next to each other). Step 2 also shows how to extract the secondary photon production spectrum in avalanche, which is responsible both for the internal and external crosstalk. To perform an internal crosstalk simulation the users can combine steps 2 and 3.

Step 2 shows how to setup a dipole source that models the secondary emission, while step 3 shows how to setup an absorption box to collect the crosstalk in neighboring SPADs. By combining a dipole source in the triggered SPAD and absorption boxes in neighboring SPADs it is in principle possible to simulate internal crosstalk as a probability to trigger a neighboring SPAD if the original SPAD is triggered. To obtain the probability of triggering the neighboring SPADs from the absorbed power, the users can use the generation rate result instead, multiply it by the avalanche triggering probability and integrate over the absorption box volume.

If the pixel size is large (more than a few microns) a full wave electromagnetic simulation with 3D FDTD may require large computational resources, especially since multiple SPADs need to be included in the simulation domain for internal crosstalk simulations. If such computational resources are not available in terms of memory and number of CPUs, the users may try to do a 2D FDTD simulation along different axes, or combine FDTD with raytracing tools, such as Ansys Zemax and Speos.

Different frequencies have different absorption lengths, so the users may choose to run multiple simulations with different simulation domain thickness in Si substrate to save time: shorter frequencies will be absorbed at shallower locations and may not contribute to the direct optical crosstalk due to reflections from the bottom of the substrate. For shorter frequencies one may choose a thinner simulation domain and focus on crosstalk though upper layers (e.g. through SiO2 in this example).

## Additional Resources

Additional documentation, examples and training material

### Related Publications

- J. B. McLaughlin, G. Gallina, F. Retière, A. De St. Croix, P. Giampa, M. Mahtab, P. Margetak, L. Martin, N. Massacret, J. Monroe, M. Patel, K. Raymond, J. Roiseux, L. Xie and G. Zhang, "Characterisation of SiPM Photon Emission in the Dark," Sensors, vol. 21, no. 17, pp. 1-15, 2021.
- G. Gallina, F. Retière, P. Giampa, J. Kroeger, P. Margetak, S. B. Mamahit, A. De St. Croix, F. Edaltafar, L. Martin, N. Massacret, M. Ward and G. Zhang, "Characterization of SiPM Avalanche Triggering Probabilities," IEEE Transactions on Electron Devices, vol. 66, no. 10, pp. 4228-4234, 2018.

### See Also

- Far field projections in FDTD overview – Lumerical Support
- FDTD vs stackdipole - halfspace emission in a multilayer stack – Lumerical Support
- Integrating power in far field projections – Lumerical Support
- Incoherent unpolarized dipole – Lumerical Support
- Avalanche photodetector – Lumerical Support
- Vertical photodetector – Lumerical Support

## Appendix

Additional background information and theory

### Finding the Optimal Thickness

The STACK simulation is chosen to initially evaluate the thickness due to its short simulation time. The intermediate steps conducted in Step 1 are shown here. The peak positions calculated in step 1 are derived from horizontal slices of results shown here.

This plot shows the change in interference pattern (transmission probability) with SiO2 thickness as simulated with STACK.