In this example, the angular response of a CMOS image sensor is characterized through optical simulations using the FDTD solver and electrical simulations using the CHARGE solver. Key results from the simulations include the spatial field profiles, transmission and optical efficiency vs. angle, quantum efficiency vs. angle. The effect of microlens shift is also considered.
Overview
Understand the simulation workflow and key results
Characterization of CMOS images sensors generally requires both optical and electrical simulations to account for the absorption, scattering, and diffraction from sub-wavelength features as well as the electrical transport of generated charge. In this example, optical simulations provide information about field profile, transmission, optical efficiency. The effects of injection angle and the microlens shift are also considered. Steps 1-3 demonstrate a few example tasks with increasing complexity (single simulation, angle/polarization sweep and angle/polarization/microlens position sweep).The generated charge data from the optical simulations (step 2) is combined with the weighting function from an independent electrical simulation (step 4) for further calculation of the quantum efficiency and the crosstalk in terms of the injection angle (step 5).
Note that the definition of the “pixel” can differ depending on the application areas. The optical simulations in this example contain a periodic array of Red/Green/Blue/Green unit cells. Throughout this example, we will refer to each of the R/G/B/G regions as a “pixel”, meaning there are 4 pixels in a unit cell, as shown in the figure below.
As this example requires many sweeps, we limited ourselves to single frequency simulations to reduce the overall simulation time. But the approaches in the optical simulations are applicable to broadband simulations.
Step 1: Initial simulation
Obtain the field profile, transmission and optical efficiency of each pixel when the sensor is illuminated by a planewave at a fixed angle. The main purpose of this step is to ensure the simulation is set up correctly and to allow the user to manually explore the results, before running the full angular response parameter sweep in the latter steps.
Step 2: Angular response
Calculate the optical efficiency and the electron-hole pair generation rates as a function of injection angle. In this example, the generation rate results are averaged in the y-direction and saved in a 2D format so that it can be used in step 5 to calculate the quantum efficiency of the device.
Step 3: Effect of microlens shift
Obtain a 2D map of the optical efficiency as a function of the injection angle and microlens shift. This sweep demonstrates an interesting way to optimize the optical performance of the device. In the interest of keeping the example short, we don’t record the generation rate data as in this step (although we could).
Step 4: Weighting function
Run the CHARGE solver to get the impulse response (Green’s function) of the system to an electron-hole pair at arbitrary positions in the substrate. From this, we calculate a spatially-varying weighting function which represents the probability of an electron-hole pair generated at any point in space to be collected by the contact for a specific pixel (green in this example). 2D CHARGE simulations are used to reduce the simulation time, but it is possible to extend the simulation methodology to 3D.
Step 5: Quantum efficiency and crosstalk
The weighting function (Step 4) is multiplied by the generation rate data (Step 2) and integrated to yield the internal quantum efficiency (IQE) and crosstalk. This approach based on the Green’s function is very efficient in calculating the IQE for arbitrary optical generation rate profiles since it requires, on the electrical side, only a single simulation for the weighting function calculation. For the sake of keeping the simulation time short, we ran a 2D CHARGE simulation to obtain a 2D weighting function and used the generation rate data that was averaged in the y-direction.
Run and results
Instructions for running the model and discussion of key results
Step 1: Initial simulation
- Open the simulation file (CMOS_image_sensor_angular_response.fsp)
- Run the script file (CMOS_image_sensor_angular_response_initial.lsf) to run the simulation and visualize some of the representative results shown below.
Field profile
The “field_XZ” and “field_YZ” frequency monitors record the fields at the cross sections of the red-green pixels and the green-blue pixels, respectively. As the source is currently set to emit at 550 nm (green), a high transmission is observed at the green pixel due to the different wavelength-selective filters on different regions.
The movies of propagating fields at the same positions as the frequency monitors can be found in the folder where the simulation file is saved. They clearly show that the injected light is selectively transmitted through the green filter and eventually dissipated by absorption in the underlying silicon layer.
Transmission
The “pixel_transmission” analysis group records the normal component of the Poynting vector, \(P_{z}(x,y)\), on the top surface of the Si layer. To calculate the power absorbed in each pixel, (optical efficiency), we can choose to integrate Pz only over the depletion region of the pixel. The easiest way to integrate \(P_{z}\) over an arbitrary region is to use a spatial filter in the shape of the depletion region and multiply it to the \(P_{z}\) . The spatial filter is optional and can be disabled in the Analysis->Script tab. The following figures show the unfiltered \(P_{z}\) , the depletion regions and the \(P_{z}\) in the depletion regions. The shape of each depletion region is currently set to a 1x1um square with a rounded corner. See integrating the poynting vector for more information.
Optical efficiency
Optical efficiency is defined as the fraction of the power incident in the pixel that is absorbed in the depletion region of the pixel:
$$\text{Optical efficiency (OE)} =\frac{\text{Absorbed power}}{\text{Source power}}$$
By integrating the \(P_{z}\) over the entire surface of the Si layer and normalizing it by the injected power, we find that about 38% of the power is transmitted into the Si layer. The combined efficiency of the two green pixels makes about 33 % while the efficiency for the red and blue pixel is about 0.5 % each.
Power into Si layer: 0.372
Power through red pixel: 0.00459
Power through 2 green pixels: 0.328
Power through blue pixel: 0.00459
The simulation file comes with a pre-run "convergence" sweep object, which records the optical efficiencies vs. mesh accuracy. The figure below shows that there relatively small differences in the results for different mesh accuracies, with the mesh accuracy of "1" giving a reasonably close result to the one for mesh accuracy "6". Using a coarse mesh for initial simulations is strongly recommended due to the large time savings it provides. The default mesh accuracy in this example is "2". However, a mesh accuracy of "1" is used in Step 3 due to the large number of simulations required. Like all simulations, thorough convergence testing is required to ensure the results are accurate.
Step 2: Angular response
- Open the simulation file (CMOS_image_sensor_angular_response.fsp)
- Open the script file (CMOS_image_sensor_angular_response_sweep_angle.lsf) and set the value of the "Analysis_only" parameter - "1" to visualize the pre-run sweep results or "0" to run the sweep and visualize the results.
- Run the script file.
The “angle sweep” consists of 14 sweep points – 7 for the source angle and 2 for the polarization. To reduce the simulation time, unnecessary simulation objects like movie monitors and index monitors were disabled.
Optical efficiency
The optical efficiency vs. source angle for each pixel is shown below. Note that we averaged the efficiencies for 0 and 90 (degree) polarizations to obtain the response for an unpolarized light. Note that a factor of \(cos(\theta)\) was also multiplied to correct for the varying amount of incident power per pixel as a function of angle. Even with an ideal pixel design, it would not be possible to have an optical efficiency as a function of angle that does better than \(cos(theta)\). The OE for green has maximum at normal incidence and is reduced at larger incidence angle. The angular response simulations also provide a measure of optical crosstalk – some light being absorbed in the red or blue pixels under green illumination (or vice versa).
Generation rate
The generation rate data from "CW_generation_gb" analysis group is used in the electrical simulation in step 4. Once the sweep is completed, 14 data files for the generation rate of the green/blue pixels will be created under the folder named "sweepangle". File names are automatically assigned by the model setup script based on the polarization angle and the source angle. The figure below shows the generation rate for an unpolarized light (550 nm) in the green/blue pixels. Note that the CW generation analysis group is currently set to average the generation rate in the "y" direction and produce an averaged 2d map of it, \(G_{L}(x,z)\). This is to make it compatible with the 2D CHARGE simulations in step 4. For the sake of saving the simulation time, we are using 2D CHARGE simulations in this example. However, it is possible to do full 3D simulations in CHARGE and hence use the full 3D generation rate from FDTD.
Step 3: Microlens shift
- Open the simulation file (CMOS_image_sensor_angular_response.fsp)
- Open the script file (CMOS_image_sensor_angular_response_microlens_shift.lsf) and set the value of the "Analysis_only" parameter - "1" to visualize the pre-run sweep results or "0" to run the sweep and visualize the results.
- Run the script file.
The "microlens shift" sweep consists of a total of 462 sweep points – 21 for microlens shift, 7 for angle theta and 2 for polarization angle. Considering the huge number of sweeps, we are using the mesh accuracy of "1" for this example to reduce the simulation time.
Optical efficiency
The optical efficiency in terms of the angle and the lens shift is shown below for each pixel. Note that the script applies a factor of \(cos(\theta)\) to correct for the varying amount of incident power per pixel at different injection angle. From the result for green, it can be seen that a shift of about 37 (nm/degree), marked by the dotted black line, gives the maximum optical efficiency for a given incidence angle. For example, if the light is incident dominantly at 15 degree, the lens needs to be shifted by about 555 nm for maximum efficiency.
Step 4: Weighting function
Note: The generation rate data from Step 2 is a prerequisite to step 4.
- Open the simulation file (CMOS_image_sensor_greens_function.ldev) and run the simulation
- Run the script file (CMOS_image_sensor_weighting_function.lsf)
Silicon has a negligible thermal recombination at moderate illumination strengths. With it being an indirect band gap semiconductor, its radiative recombination (generation of photons by recombination of electron-hole pairs) is also negligible. As a result, the photo-generated charges in the substrate will be mostly collected by the contacts for different pixels. In this step, we use the point generation source in CHARGE simulation to determine the weighting function, \(W(x,y,z)\) of the device. \( W(x,y,z)\) is the probability that the charge generated at that position will be collected by a specific contact. This approach is based on the Green’s function, \(G(x,y,z)\), and is equivalent to knowing the carrier density \(n,p\) at each contact in response to an impulse generation rate source located at a specific location.
The carrier densities are calculated during the charge transport simulation. To determine the full \(G(x,y,z)\), the location of the impulse source is swept over all r in the simulation domain. This operation is performed internally when the “map current collection probability” option is enabled for one or more contacts. When the simulation is complete, the weighting functions are stored as the result "W" for each carrier type at each enabled contact and are exported as data files for further analysis in step 5. The figure below shows the weighting function \(W(x,z,)\) for the green pixel shows that the collection probability for the green pixel is very high when the charge is located nearer the green contact (top left). However, it also shows that some of the charges generated in the blue pixel region (\(x>0\)) has a non-zero probability of being collected by the green contact. This suggests that there is some electrical crosstalk between the neighbouring pixels.
Step 5: IQE and crosstalk
Note: The generation rate data from Step 2 is a prerequisite to step 5.
- Open the simulation file (CMOS_image_sensor_greens_function.ldev)
- Run the script file (CMOS_image_sensor_angular_response_iqe.lsf)
In this step, we will be calculating the quantum efficiency (QE) of the green pixel and the green/blue crosstalk based on the Green’s function approach. We do not need to run any additional simulations at this stage as we already have all the required generation rate and weighting function data saved from the step 2 and 4, respectively. The definitions of the related quantities are as follows:
$$\text{Internal quantum efficiency (IQE)} = \frac{\text{Charges collected by the green pixel}}{\text{Total charges generated by absorption}}$$
$$\text{External quantum efficiency (EQE)} = \text{IQE}\times \text{OE (Optical efficiency)}$$
$$\text{Green/blue crosstalk} = \frac{\text{Charges collected by the blue pixel}}{\text{Total charges generated by absorption}}$$
Quantum efficiency and crosstalk
The script sequentially loads the 14 generation rate data saved from the angle sweep in step 2 and multiplies it with the weighting function for green. The following figures show \(G_{L}(x,z)\), \(W_{green}(x,z)\) and \(G_{L}(x,z)W_{green}(x,z)\) for the unpolarized light at normal incidence.
By integrating the \(G_{L}(x,z)W_{green}(x,z)\) and normalizing it with the total generation rate, we obtain the IQE for the green pixel. Repeating the same procedure with the weighting function for the blue pixel, \(W_{blue}(x,z)\), yields the green/blue crosstalk. The IQE has its maximum value of about 80 % and decreases at larger source angle. The trend is in agreement with the increased green/blue crosstalk at larger angle. The maximum EQE is about 26%. When interpreting this data, it’s important to remember that \(G_{L}(x,z)\) is for green light illumination.
Important model settings
Description of important objects and settings used in this model
Parametrization
All the structures in this example are constructed using the setup script in the “image sensor” structure group. However, some of the key design parameters such as pixel size and the microlens shift are set in the setup tab of the “model”, which then will update the associated parameters in the “image sensor” as well as in other simulation objects dependent on any of those parameters.
"CW generation" analysis groups
The generation rate group consists of a 3D frequency monitor. It measures the absorbed power in the Si layer, then calculates the generation rate assuming a single pair of electron-hole pair is created per single photon.
- Source propagation axis: The script in the “CW generation” group assumes the source is injecting in the y-direction for 2D and the z-direction for 3D simulations.
- Averaging of generation rate: Even though the raw data for the generation rate is obtained in 3D, we are averaging it in the y-direction and saves a 2D generation data for later use. This is because we are running a 2D simulation in CHARGE to save the simulation time in this demonstrative example.
- x/y/z spans: The spans of the generation rate objects are set by the setup script in the “model.” The z-min of the objects might need to be adjusted to capture most of the absorbed light penetrating deeper into the substrate. The z-max of the objects needs to be at least one mesh cell away from the Si surface for an accurate absorption calculation . You might need to do some convergence test to decide on the appropriate z-span.
- Filenames: The name of the generation rate file name is automatically set by the “polarization angle” of the source and the “count” in the setup variables tab of the “model.” The “count” is paired with the parameter “theta” in the “sweep angle” object. At each sweep point of “sweep angle”, the “count” value in “model” is updated and consequently the “export filename” in the CW generation group.
Unpolarized light
To obtain the transmission and generation rate results for unpolarized light, we need to run two simulations (one for x-polarized and another for y-polarized light) and average them. For that reason, a nested sweep named “source polarization” is included in the “angle sweep” and the “microlens shift” sweep objects.
For initial simulations, it can make sense to ignore polarization effects and choose only one polarization, S or P. Eventually, it makes sense to correctly calculate the incoherent response using both polarizations but a great deal of initial testing and optimization can be done with only one polarization.
"map current collection probability"
To calculate the current collection probability (=weighting function) for a given electrical contact, select the “map current collection probability” option in the steady-state contact configuration. In this simulation, the collection probability mapping is enabled for each of the simulation contacts located in the PD(Photodiode) n-wells. These contacts are biased to mimic the potential of a depletion region. The p++ surface diffusion and the substrate are electrically grounded.
Use a coarse mesh initially
We strongly recommend starting with a coarse mesh, using a mesh accuracy setting of 1 or 2. It is far easier to setup, test or optimize a simulation that takes a few seconds to minutes, compared to a simulation that takes several hours. Only once all other problems are resolved and you are getting good results should you try using a mesh accuracy setting of 3 or more by doing some convergence testing.
X,Y override in Si
Due to the high index of Silicon (n>4), the automatic meshing algorithm will use a very small mesh everywhere in the Si region. This small mesh will make the simulation significantly slower that it would otherwise be.
From Snell’s law, we can show that the in-plane wave vector is conserved (\(k_{x, glass} = k_{x, Si}\)).
In other words, the in-plane wavevector (kx or ky) in the Silicon can never be larger than the in-plane wavevector in the glass. Therefore, it is possible to use a coarser mesh in the x and y directions without reducing the accuracy. To force a larger mesh size in the X and Y directions, add a mesh override region over the Si layer. Set the override region to treat this area as if it has an index of 1.5 (glass) in the X and Y directions. A small mesh is required in the Z direction, so the override should not be applied in the Z direction.
This technique is only valid if there are no scattering structures within the Si. Indeed, even scattering structure on or near the surface of the Si can generate light propagating at all angles in the Si, however, the above approximation is generally very good even if there are some scattering structures on the surface of the Si. We recommend testing the convergence if there are concerns about the amount of steep angle scattering that may be generated in the Si.
Use PEC for metal objects
In most image sensor simulations, metallic objects behave like perfect metals (100% reflection, no absorption). Rather than using a detailed material model for the metal (which will require a smaller mesh), simply use the Perfect Electrical Conductor (PEC) material model. The PEC material model does not require a small mesh, which makes your simulations faster.
Use conformal meshing (Conformal variant 1)
The conformal mesh can provide much more accurate results at larger mesh sizes and make it possible to run simulation much faster. If PEC is used for the metals, it is a good idea to switch to using the "Conformal variant 1" setting which will apply the conformal mesh algorithm to the PEC as well as the dielectric materials and the Si (please see Mesh refinement options for more details). The figure on the right shows that the Optical Efficiency to the green pixel is almost unchanged when going from a mesh accuracy of 1 (lambda/dx=6) to a mesh accuracy of 6 (lambda/dx=26). When using staircase meshing, the convergence is slower and a smaller mesh size is required for the same accuracy.
Please note that the conformal mesh can generate more numerical instability, especially if many highly dispersive materials are used that require large numbers of coefficients to fit. If these instabilities do occur, they can normally be controlled by making changes to the fit settings for certain materials. Please contact Lumerical support for advice if necessary.
Structures in optical and electrical simulations
You do not need to include exactly the same structures and simulation volume in the optical and the electrical simulations. Structures such as microlens, color filters and metal shields do not affect the electrical response of the device and should be removed from the electrical simulations to avoid an unnecessarily larger simulation region. Likewise, some of the electrical parts such as the substrate contact do not need to be included in the optical simulation if they have negligible effect on the optical response of the device.
Nwell contacts
This contact should not exist in realistic design, but this is a key element to help calculate the weighting function. In simulation, it plays a dual role:
- set the electrostatic potential of the n-well equivalent to its depleted state after reset
- collect photo-generated charge carriers that are captured in the n-well
When calculating the collection probability weighting function, we assume that the electrostatic potential is unperturbed. In normal device operation, charge carriers trapped in this region will be transferred to the drain (via the TX gate) and will accumulate on the source-follower gate. An accurate determination of the potential of the depleted n-well can be determined by inducing a channel between the drain (biased to VDD) and the n-well with the TX gate turned ON, which can be simulated at steady state. If user would like a quick way to test the simulation behavior, they can first run a simulation without the nwell contact to test it electrostatic potential. And then add a small amount, eg, 0.5V, to this simulated voltage to represent the bias at the nwell. This contact should be small in size and placed at the minimum of the potential in the n-well (small enough that the potential is close to uniform over the equivalent of the contact’s surface). This usually corresponds to the peak location of the doping region.
Sub contact
This is a simulation contact. In realistic design, the silicon substrate should be grounded and therefore we need to assign a contact to the substrate. In the Charge solver, contacts have to be assigned to a metal type material. Aluminum is chosen arbitrarily, and the electrical contact boundary condition “force Ohmic” is enabled. The contact should be placed sufficiently far from all depletion regions so that it doesn’t change the potential there (i.e. equilibrium condition should be established before contact surface), but generally must be placed deeper in the substrate to provide adequate depth for the absorption of light (check the optical generation rate profile).
PD contacts
This simulation contact is disabled by default. If necessary, the PD objects and contacts can be enabled and used as a troubleshooting step. There are some p++ regions without a reference voltage near the surface of the silicon. The solver may become numerically unstable if the region is floating. This contact can provide a reference voltage to the region for simulation purpose. This contact should be placed on top of the p++ region, the size of this contact should not have major effects to the simulation results.
Transfer gate and drain
Those are not simulated in this example. The current example is simplified and a focus of it is to show the simulation methodologies and how the Green's function approach can be used to efficiently calculate quantities like IQE, EQE, crosstalk, etc. In a more sophisticated simulation setup, experienced users may want to include these elements in their simulations, but this will certainly increase the complexity of the simulation.
Isolation trenches
These are common features to reduce electrical leakage current, formed by etching a Shallow Trench Isolation (STI) in the substrate and filling it with oxide. Their presence and location will depend on the specific design and fabrication process.
Updating the model with your parameters
Instructions for updating the model based on your device parameters
Customizing structures
The “image sensor” structure group contains only some representative parts of typical CMOS image sensors. When modifying structures, it is recommended that you retain the “user properties” of the “image sensor” object and make sure the link between the key design parameters (pixel size, microlens shift, etc.) and the associated structures are not broken. Otherwise, you might need to write all the relevant scripts from scratch.
Importing AFM data for the lens
The “image sensor” structure group has an option for “use_AFM_data” to demonstrate how AFM data could be used to define parts of the sensor. Set its value to “1” to import the AFM data for the microlens surface and “0” to define the surface based on the polynomial and conic formulations. The script assumes you have your surface data ( savedata - Script command ) saved in the same directory as your .fsp file and uses the importsurface2 command to create the surface. You can also import a surface data in .txt format using the importsurface command. The “cmos_image_sensor_angular_response_microlens_AFM_make.lsf” can be used to generate an exemplary surface data. Please note that using AFM data will make your fsp files much larger, and you will likely want to keep a low setting for the “rendering detail” to avoid slower display of the structures.
Source wavelength
The wavelength of the source is currently set at 550nm (green). If you want to see the response of the green/blue pixels to blue light, all you need to do is just change the source wavelength to blue. However, to obtain the response of the red/green pixels to red light, you need some additional modification in the simulation settings:
- Enable the the associated generation rate group, “CW_generation_rg”
- Add the “::model::CW_generation_rg::Igen” to the “Result” of the “sweep angle” object
- Update the filenames for the weighting functions in the script files for step 4 and 5 with correct pixel names.
Source intensity
The default intensity of the light source is set to 1W/m^2 in the “model.” Its value might need to be updated to account for the actual intensity of the light source in consideration.
Spatial filter for depletion region
The depletion region is where the generated charges are swept by the strong electric field and contribute the photocurrent. The shape of the actual depletion region can differ depending on the contact design and doping profiles. In the FDTD simulation, we assumed that light absorbed within the depletion region may contribute photocurrent and therefore used a spatial filter in the analysis script of the “pixel transmission” to mimic the depletion region when calculating the transmission for each pixel. It is currently set to a 1x1 \(({\mu m}^2) \)with rounded corners, but you might need to modify it to better match the depletion region in your device. This is an optional setting.
Taking the model further
Information and tips for users that want to further customize the model
Broadband simulations
While the optical simulations in step 1-3 were run at single frequency, you can also use the same file for broadband (multi-frequency) simulations. However, there are some additional points to be considered when running broadband simulations. Please visit here for detailed information. Note that the generation rate analysis group returns results that are averaged over wavelength when broadband source is used.
3D CHARGE simulations
While the electrical part of this example was based on 2D CHARGE simulations, the Green’s function approach works for the 3D CHARGE simulations as well. Please note that extending the current example to 3D CHARGE simulations requires extensive changes to the generation rate script, the CHARGE simulation file and related analysis scripts.
Point spread function (PSF)
Together with OE, PSF is also a frequently used metric to characterize the optical properties of CMOS image sensor. Broadly speaking, it is a measure of spatial crosstalk — i.e., how much light is detected in neighboring pixels when a specific pixel is fully illuminated. Please visit here for further information and example files.
Response for arbitrary illumination
The approach in the above PSF example uses an array of “thin lens” Gaussian beams to mimic a source uniformly illuminating a specific pixel. This requires the simulation to have pml boundaries in all directions and a large simulation region. Additionally, you need to run separate simulations to obtain the results for different objective lenses. Fortunately, there is a much faster and more efficient approach available, which use an incoherent sum of results for planewave illumination to reconstruct the response of any arbitrary illumination (including arbitrary objective lens). For further information about this approach based on the planewave, please see [2]
Additional resources
Additional documentation, examples and training material
Related publications
- F. Hirigoyen, A. Crocherie, J. M. Vaillant, and Y. Cazaux, “FDTD-based optical simulations methodology for CMOS image sensors pixels architecture and process optimization” Proc. SPIE 6816, 681609 (2008)
- J. Vaillant, A. Crocherie, F. Hirigoyen, A. Cadien, and J. Pond, "Uniform illumination and rigorous electromagnetic simulations applied to CMOS image sensors," Opt. Express 15, 5494-5503 (2007)
- Crocherie et al., “Three-dimensional broadband FDTD optical simulations of CMOS image sensor”, Optical Design and Engineering III, Proc. of SPIE, 7100, 71002J (2008)
- Wang, Xinyang, "Noise in Sub-Micron CMOS Image Sensors", Ph.D. Thesis, Delft University of Technology
- W. Gazeley and D. McGrath, "Quantum Efficiency Simulation Using Transport Equations," International Image Sensor Workshop, R06 (2011)
See also
- Optical simulation methodology
- Electrical simulation methodology
- Point spread function
- Green’s function IQE method