Staring Array Module (rystare)

Overview

This module provides a high level model for CCD and CMOS staring array signal chain modelling. The model accepts an input image in photon rate irradiance units and then proceeds to calculate the various noise components and signal components along the signal flow chain.

The code in this module serves as an example of implementation of a high-level CCD/CMOS photosensor signal chain model. The model is described in the article ‘High-level numerical simulations of noise in solid-state photosensors: review and tutorial’ by Mikhail Konnik and James Welsh, arXiv:1412.4031v1 [astro-ph.IM]. The code was originally written in Matlab and used for the Adaptive Optics simulations and study of noise propagation in wavefront sensors, but can be used for many other applications involving CCD/CMOS photosensors. The original files are available at:

The original Matlab code was ported to Python and extended in a number of ways. The core of the model remains the original Konnik model as implemented in the Matlab code. The Python code was validated against results obtained with the Matlab code, up to a point and then substantially reworked and refactored. During the refactoring due diligence was applied with regression testing, checking the new results against the previous results. The results were communicated and confirmed with Konnik. A number of corrections were also made to Konnik’s original code, some with his cooperation and some without his involvement.

The documentation in the code was copied from Konnik’s Matlab code, so he deserves the credit for the very detailed documentation. His documentation was extracted from the paper quoted above.

The sample Python code (derived from Konnik’s code) in the repository models two different cases

  • a simple model: which is completely linear (no non-linearities), where all noise are Gaussian, and without source follower noise,
  • an advanced model: which has V/V and V/e non-linearities, Wald or lognormal noise, source follower and sense node noise sources and even ADC non-linearities.

The Python code supports enabling/disabling of key components by using flags.

In the documentation for the Matlab code Konnik expressed the hope “that this model will be useful for somebody, or at least save someone’s time. The model can be (and should be) criticized.” Indeed it has, thanks Mikhail! Konnik quotes George E. P. Box, the famous statistician, and who said that “essentially, all models are wrong, but some are useful”.

Signal Flow

The process from incident photons to the digital numbers appearing in the image is outlined in the picture below. The input image to rystare must be provided in photon rate irradiance units [q/(s.m2)], with photon noise already present in the image. The count of photons captured in the detector is determined from the irradiance by accounting for the detector area and integration time. Then, the code models the process of conversion from photons to electrons and subsequently to signal voltage. Various noise sources are modelled to derive at a realistic image model. Finally, the ADC converts the voltage signal into digital numbers. The whole process is depicted in the figure below.

camerascheme-horiz.png

Many noise sources contribute to the resulting noise image that is produced by the sensor. Noise sources can be broadly classified as either fixed-pattern (time-invariant) or temporal (time-variant) noise. Fixed-pattern noise refers to any spatial pattern that does not change significantly from frame to frame. Temporal noise, on the other hand, changes from one frame to the next. All these noise sources are modelled in the code. For more details see Konnik’s original paper or the docstrings present in the code.

Changes to Matlab code

  1. Renamed many, if not all, variables to be more descriptive.
  2. Created a number of new functions by splitting up the Matlab functions for increased modularity.
  3. Store (almost) all input and output variables in an HDF5 file for full record keeping.
  4. Precalculate the image data input as HDF5 files with linear detector parameters embedded in the file. This was done to support future image size calculations. The idea is to embed the target frequency in the data file to relate observed performance with the frequency on the focal plane.
  5. Moved sourcefollower calcs out from under dark signal flag. sourcefollower noise is now always calculated irrespective of whether dark noise is selected or not.
  6. Input image now photon rate irradiance q/(m2.s), image should already include photon noise in input. Removed from ccd library: irradiance from radiant to photon units, adding photon shot noise. This functionality has been added to the image generation code.
  7. Both CCD and CMOS now have fill factors, the user can set CCD fill factor differently from CMOS fill factor. The fill factor value is used as-in in the rest of the code, without checking for CCD or CMOS. This is done because CCD fill factor is 1.0 for full frame sensors but can be less than 1.0 for other types of CCD.
  8. Now uses SciPy’s CODATA constants where these are available.
  9. Put all of the code into a single file rystare.py in the pyradi repository.
  10. Minor changes to Konnik’s excellent documentation to be Sphinx compatible. Documentation is now generated as part of the pyradi documentation.
  11. The original model did not implement detector shot noise. The absence of shot noise has a huge effect for low-electron-count signals. The detector shot noise is now correctly modelled here.

Example Code

The two examples provided by Konnik are merged into a single code, with flags to select between the two options. The code is found at the end of the module file in the __main__ part of the module file. Set doTest = ‘Simple’ or doTest = ‘Advanced’ depending on which model. Either example will run the photosensor function (all functions are thoroughly documented in the Python code, thanks Mikhail!).

The two prepared image files are both 256x256 in size. New images can be generated following the example shown in the __main__ part of the rystare.py module file (use the function create_HDF5_image as a starting point to develop your own).

The easiest way to run the code is to open a command window in the installation directory and run the run_example function in the module code. This will load the module and execute the example code function. Running the example code function will create files with names similar to PSOutput.hdf5 and PSOutput.txt. To run the example, create a python file with the following contents and run it at the command line prompt:

import pyradi.rystare as rystare
rystare.run_example('Advanced','Output', doPlots=True, doHisto=True, doImages=True)
rystare.run_example('Simple','Output', doPlots=True, doHisto=True, doImages=True)

By setting all the flags to True the example code will print a number of images to file. Plotting the results to file takes a while. Execution is much faster with all flags set to False.

Study the text file using a normal text editor and study the HDF5 file by using the viewer available from https://www.hdfgroup.org/products/java/hdfview/.

An IPython notebook demonstrating this model in quite some detail is available here https://github.com/NelisW/ComputationalRadiometry/blob/master/09b-StaringArrayDetectors.ipynb

The full code for the example run is included in rystare.py:

#prepare so long for Python 3
from __future__ import division
from __future__ import print_function
from __future__ import unicode_literals

import numpy as np
import re
import os.path
from matplotlib import cm as mcm
import matplotlib.mlab as mlab

import pyradi.rystare as rystare
import pyradi.ryplot as ryplot
import pyradi.ryfiles as ryfiles
import pyradi.ryutils as ryutils

     ################################################################
    def run_example(doTest='Advanced', outfilename='Output', pathtoimage=None,
        doPlots=False, doHisto=False, doImages=False):
        """This code provides examples of use of the pyradi.rystare model for
        a CMOS/CCD photosensor.

        Two models are provided 'simple' and 'advanced'

        doTest can be 'Simple' or 'Advanced'


        Args:
            | doTest (string):  which example to run 'Simple', or 'Advanced'
            | outfilename (string):  filename for output files
            | pathtoimage (string):  fully qualified path to where the image is located
            | doPlots (boolean):  flag to control the creation of false colour image plots with colour bars
            | doHisto (boolean):  flag to control the creation of image histogram plots
            | doImages (boolean):  flag to control the creation of monochrome image plots

        Returns:
            | hdffilename (string): output HDF filename

        Raises:
            | No exception is raised.

        Author: Mikhail V. Konnik, revised/ported by CJ Willers

        Original source: http://arxiv.org/pdf/1412.4031.pdf
        """
        import os.path
        from matplotlib import cm as mcm
        import matplotlib.mlab as mlab

        import pyradi.ryfiles as ryfiles
        import pyradi.ryutils as ryutils


        if doTest in ['Simple']:
            prefix = 'PS'
        elif  doTest in ['Advanced']:
            prefix = 'PA'
        else:
            exit('Undefined test')

        [m, cm, mm, mum, nm, rad, mrad] = define_metrics()

        #open the file to create data structure and store the results, remove if exists

        hdffilename = '{}{}.hdf5'.format(prefix, outfilename)
        if os.path.isfile(hdffilename):
            os.remove(hdffilename)
        strh5 = ryfiles.open_HDF(hdffilename)

        # Light Noise parameters
        strh5['rystare/photonshotnoise/activate'] = True #photon shot noise.

        #sensor parameters
        strh5['rystare/sensortype'] = 'CCD' # CCD / CMOS must be in capitals

        # full-frame CCD sensors has 100% fil factor (Janesick: 'Scientific Charge-Coupled Devices')
        if strh5['rystare/sensortype'].value in ['CMOS']:
            strh5['rystare/photondetector/geometry/fillfactor'] = 0.5 # Pixel Fill Factor for CMOS photo sensors.
        else:
            strh5['rystare/photondetector/geometry/fillfactor'] = 1.0 # Pixel Fill Factor for full-frame CCD photo sensors.

        strh5['rystare/photondetector/integrationtime'] = 0.01 # Exposure/Integration time, [sec].
        strh5['rystare/photondetector/externalquantumeff'] = 0.8  # external quantum efficiency, fraction not reflected.
        strh5['rystare/photondetector/quantumyield'] = 1. # number of electrons absorbed per one photon into material bulk


        # photo response non-uniformity noise (PRNU), or also called light Fixed Pattern Noise (light FPN)
        strh5['rystare/photondetector/lightPRNU/activate'] = True
        strh5['rystare/photondetector/lightPRNU/seed'] = 362436069
        strh5['rystare/photondetector/lightPRNU/model'] = 'Janesick-Gaussian'
        strh5['rystare/photondetector/lightPRNU/sigma'] = 0.01 # sigma [about 1\% for CCD and up to 5% for CMOS]

        # detector material bandgap properties
        strh5['rystare/photondetector/varshni/Egap0'] = 1.166  #bandgap energy for 0 degrees of K. [For Silicon, eV]
        strh5['rystare/photondetector/varshni/varA'] = 5.5e-04 #Si material parameter, [eV/K].
        strh5['rystare/photondetector/varshni/varB'] = 636. #Si material parameter, [K].

        # Dark Current Noise parameters
        strh5['rystare/photondetector/darkcurrent/activate'] = True
        strh5['rystare/photondetector/operatingtemperature'] = 300. # operating temperature, [K]
        strh5['rystare/photondetector/darkcurrent/ca'] = 4.31e5 # for density in m2
        strh5['rystare/photondetector/darkcurrent/ed'] = 2.
        strh5['rystare/photondetector/darkcurrent/densityAm2'] = 1. # dark current figure of merit, [nA/cm2].  For very poor sensors, add DFM
        #  Increasing the DFM more than 10 results to (with the same exposure time of 10^-6):
        #  Hence the DFM increases the standard deviation and does not affect the mean value.

        # dark current shot noise
        strh5['rystare/photondetector/darkcurrent/shotnoise/activate'] = True
        strh5['rystare/photondetector/darkcurrent/shotnoise/seed'] = 6214069
        strh5['rystare/photondetector/darkcurrent/shotnoise/model'] = 'Gaussian'

        #dark current Fixed Pattern Noise
        strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/activate'] = True
        # Janesick's book: dark current FPN quality factor is typically between 10\% and 40\% for CCD and CMOS sensors
        strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/seed'] = 362436128

        if doTest in ['Simple']:
            strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/model'] = 'Janesick-Gaussian'
            strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/sigma'] =  0.3 #0.3-0.4 sigma for dark current signal (Janesick's book)
        elif  doTest in ['Advanced']:
            strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/model'] = 'LogNormal' #suitable for long exposures
            strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/sigma'] = 0.4 # lognorm_sigma.
        else:
            pass

        # #alternative model
        # strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/model']  = 'Wald'
        # strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/sigma']  = 2.0 #small parameters (w<1) produces extremely narrow distribution, large parameters (w>10) produces distribution with large tail.

        # #alternative model
        # strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/model']  = 'AR-ElGamal'
        # strh5['rystare/photondetector/darkcurrent/fixedPatternNoise/filter_params']  = [1., 0.5] # see matlab filter or scipy lfilter functions for details

        #sense node charge to voltage
        strh5['rystare/sensenode/gain'] = 5e-6 # Sense node gain, A_SN [V/e]
        strh5['rystare/sensenode/vrefreset'] = 3.1 # Reference voltage to reset the sense node. [V] typically 3-10 V.
        strh5['rystare/sensenode/vsnmin'] = 0.5 # Minimum voltage on sense node, max well charge [V] typically < 1 V.
        strh5['rystare/sensenode/gainresponse/type'] = 'linear'
        strh5['rystare/sensenode/gainresponse/k1'] = 1.090900000e-14 # nonlinear capacitance is given by C =  k1/V
        if strh5['rystare/sensenode/gainresponse/type'] in ['nonlinear']:
            strh5['rystare/sensenode/fullwellelectronselection/fullwellelectrons'] = \
                -(strh5['rystare/sensenode/gainresponse/k1'].value/const.e) * \
                np.log(strh5['rystare/sensenode/vsnmin'].value/strh5['rystare/sensenode/vrefreset'].value)
        else:
            strh5['rystare/sensenode/fullwellelectronselection/fullwellelectrons'] = 2e4 # full well of the pixel (how many electrons can be stored in one pixel), [e]

        strh5['rystare/sensenode/resetnoise/activate'] = True
        strh5['rystare/sensenode/resetnoise/factor'] = 0.8 # the compensation factor of the Sense Node Reset Noise:
                                               # 1 - no compensation from CDS for Sense node reset noise.
                                               # 0 - fully compensated SN reset noise by CDS.
        strh5['rystare/sensenode/resetnoise/seed'] = 2154069
        strh5['rystare/sensenode/resetnoise/model'] = 'Gaussian'


        #source follower
        strh5['rystare/sourcefollower/gain'] = 1. # Source follower gain, [V/V], lower means amplify the noise.


        #source follower
        strh5['rystare/sourcefollower/nonlinearity/activate'] = True # VV non-linearity
        strh5['rystare/sourcefollower/nonlinearity/ratio'] = 1.05 # > 1 for lower signal, < 1 for higher signal
        strh5['rystare/sourcefollower/noise/flickerCornerHz'] = 1e6 #flicker noise corner frequency $f_c$ in [Hz], where power spectrum of white and flicker noise are equal [Hz].
        strh5['rystare/sourcefollower/noise/whitenoisedensity'] = 15e-9 #thermal white noise [\f$V/Hz^{1/2}\f$, typically \f$15 nV/Hz^{1/2}\f$ ]
        strh5['rystare/sourcefollower/noise/deltaindmodulation'] = 1e-8 #[A] source follower current modulation induced by RTS [CMOS ONLY]
        strh5['rystare/sourcefollower/dataclockspeed'] = 20e6 #MHz data rate clocking speed.
        strh5['rystare/sourcefollower/freqsamplingdelta'] = 10000. #sampling spacing for the frequencies (e.g., sample every 10kHz);
        if doTest in ['Simple']:
            strh5['rystare/sourcefollower/noise/activate'] = False
        elif  doTest in ['Advanced']:
            strh5['rystare/sourcefollower/noise/activate'] = True

        #dark current Offset Fixed Pattern Noise
        strh5['rystare/sourcefollower/fpoffset/activate'] = True
        strh5['rystare/sourcefollower/fpoffset/model'] = 'Janesick-Gaussian'
        strh5['rystare/sourcefollower/fpoffset/sigma'] = 0.0005 # percentage of (V_REF - V_SN)
        strh5['rystare/sourcefollower/fpoffset/seed'] = 362436042


        # Correlated Double Sampling (CDS)
        if doTest in ['Simple']:
            strh5['rystare/sourcefollower/CDS/sampletosamplingtime'] = 0 #not used
        elif  doTest in ['Advanced']:
            strh5['rystare/sourcefollower/CDS/sampletosamplingtime'] = 1e-6 #CDS sample-to-sampling time [sec].
        else:
            pass
        strh5['rystare/sourcefollower/CDS/gain'] = 1. # CDS gain, [V/V], lower means amplify the noise.

        # Analogue-to-Digital Converter (ADC)
        strh5['rystare/ADC/num-bits'] = 12. # noise is more apparent on high Bits
        strh5['rystare/ADC/offset'] = 0. # Offset of the ADC, in DN
        strh5['rystare/ADC/nonlinearity/activate'] = False
        strh5['rystare/ADC/nonlinearity/ratio'] = 1.1

        #Sensor noises and signal visualisation
        strh5['rystare/flag/plots/doPlots'] = False
        strh5['rystare/flag/plots/plotLogs'] = False



        #For testing and measurements only:
        strh5['rystare/darkframe'] = False # True if no signal, only dark

        #=============================================================================

        if strh5['rystare/darkframe'].value:  # we have zero light illumination
            imagehdffilename = 'data/image-Zero-256-256.hdf5'
        else:   # load an image, nonzero illumination
            imagehdffilename = 'data/image-Disk-256-256.hdf5'
            # imagehdffilename = 'data/image-Uniform-256-256.hdf5'

        if pathtoimage is None:
            pathtoimage = os.path.join(os.path.dirname(ryprob.__file__), imagehdffilename)

        imghd5 = ryfiles.open_HDF(pathtoimage)

        #images must be in photon rate irradiance units q/(m2.s)
        strh5['rystare/equivalentSignal'] = imghd5['image/equivalentSignal'].value
        strh5['rystare/signal/photonRateIrradianceNoNoise'] = imghd5['image/PhotonRateRadianceNoNoise'].value
        strh5['rystare/signal/photonRateIrradiance'] = imghd5['image/PhotonRateRadiance'].value
        strh5['rystare/pixelPitch'] = imghd5['image/pixelPitch'].value
        strh5['rystare/imageName'] = imghd5['image/imageName'].value
        strh5['rystare/imageFilename'] = imghd5['image/imageFilename'].value
        strh5['rystare/imageSizePixels'] = imghd5['image/imageSizePixels'].value
        strh5['rystare/wavelength'] = imghd5['image/wavelength'].value
        strh5['rystare/imageSizeRows'] = imghd5['image/imageSizeRows'].value
        strh5['rystare/imageSizeCols'] = imghd5['image/imageSizeCols'].value
        strh5['rystare/imageSizeDiagonal'] = imghd5['image/imageSizeDiagonal'].value
        strh5['rystare/equivalentSignalUnit'] = imghd5['image/equivalentSignalUnit'].value
        strh5['rystare/LinUnits'] = imghd5['image/LinUnits'].value



        #calculate the noise and final images
        strh5 = photosensor(strh5) # here the Photon-to-electron conversion occurred.

        with open('{}{}.txt'.format(prefix,outfilename), 'wt') as fo:
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('SignalPhotonRateIrradiance',np.mean(strh5['rystare/signal/photonRateIrradiance'].value), np.var(strh5['rystare/signal/photonRateIrradiance'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('signalphotonRateIrradianceNU',np.mean(strh5['rystare/signal/photonRateIrradianceNU'].value), np.var(strh5['rystare/signal/photonRateIrradianceNU'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('signalelectronRateIrradiance',np.mean(strh5['rystare/signal/electronRateIrradiance'].value), np.var(strh5['rystare/signal/electronRateIrradiance'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('SignalelectronRate',np.mean(strh5['rystare/signal/electronRate'].value), np.var(strh5['rystare/signal/electronRate'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('signallightelectronsNoShotNoise',np.mean(strh5['rystare/signal/lightelectronsnoshotnoise'].value), np.var(strh5['rystare/signal/lightelectronsnoshotnoise'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('signallightelectrons',np.mean(strh5['rystare/signal/lightelectrons'].value), np.var(strh5['rystare/signal/lightelectrons'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('signalDark',np.mean(strh5['rystare/signal/darkcurrentelectrons'].value), np.var(strh5['rystare/signal/darkcurrentelectrons'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('source_follower_noise',np.mean(strh5['rystare/sourcefollower/source_follower_noise'].value), np.var(strh5['rystare/sourcefollower/source_follower_noise'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('SignalElectrons',np.mean(strh5['rystare/signal/electrons'].value), np.var(strh5['rystare/signal/electrons'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('voltagebeforeSF',np.mean(strh5['rystare/signal/voltagebeforeSF'].value), np.var(strh5['rystare/signal/voltagebeforeSF'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('voltagebeforecds',np.mean(strh5['rystare/signal/voltagebeforecds'].value), np.var(strh5['rystare/signal/voltagebeforecds'].value)))
            # fo.write('{:26}, {:.5e}, {:.5e}\n'.format('voltagebeforeadc',np.mean(strh5['rystare/signal/voltagebeforeadc'].value), np.var(strh5['rystare/signal/voltage'].value)))
            # fo.write('{:26}, {:.5e}, {:.5e}\n'.format('SignalVoltage',np.mean(strh5['rystare/signal/voltage'].value), np.var(strh5['rystare/signal/voltagebeforeadc'].value)))
            fo.write('{:26}, {:.5e}, {:.5e}\n'.format('SignalDN',np.mean(strh5['rystare/signal/DN'].value), np.var(strh5['rystare/signal/DN'].value)))

        lstimgs = ['rystare/signal/photonRateIrradianceNoNoise', 'rystare/quantumEfficiency',
             'rystare/signal/photonRateIrradiance','rystare/photondetector/lightPRNU/value',
             'rystare/signal/photonRateIrradianceNU','rystare/signal/electronRateIrradiance',
             'rystare/signal/electronRate', 'rystare/signal/lightelectronsnoshotnoise','rystare/signal/lightelectrons',
             'rystare/darkcurrentelectronsnonoise', 'rystare/signal/darkcurrentelectrons',
             'rystare/photondetector/darkcurrent/fixedPatternNoise/value',
             'rystare/signal/darkcurrentelectrons',
             'rystare/signal/electrons','rystare/signal/electronsWell',
             'rystare/signal/sensenodevoltageLinear','rystare/noise/sn_reset/resetnoise',
             'rystare/noise/sn_reset/vrefresetpluskTC','rystare/sensenode/vrefreset',
             'rystare/signal/voltage','rystare/sourcefollower/gainA','rystare/signal/voltageAfterSF',
             'rystare/sourcefollower/source_follower_noise','rystare/signal/voltageAfterSFnoise',
             'rystare/sourcefollower/fpoffset/value','rystare/signal/voltagebeforecds',
             'rystare/signal/voltageaftercds','rystare/ADC/gainILE','rystare/ADC/gain','rystare/signal/DN']

        if doPlots:
            ryfiles.plotHDF5Images(strh5, prefix=prefix, colormap=mcm.jet,  lstimgs=lstimgs,
                logscale=strh5['rystare/flag/plots/plotLogs'].value, debug=False)

        if doHisto:
            ryfiles.plotHDF5Histograms(strh5, prefix, bins=100, lstimgs=lstimgs)

        if doImages:
            ryfiles.plotHDF5Bitmaps(strh5, prefix, pformat='png', lstimgs=lstimgs,debug=False)

        strh5.flush()
        strh5.close()

        return hdffilename

HDF5 File

The Python implementation of the model uses an HDF5 file to capture the input and output data for record keeping or subsequent analysis. HDF5 files provide for hierarchical data structures and easy read/save to disk. See the file hdf5-as-data-format.md ([hdf5asdataformatalt]) in the pyradi root directory for more detail.

Input images are written to and read from HDF5 files as well. These files store the image as well as the images’ dimensional scaling in the focal plane. The intent is to later create test targets with specific spatial frequencies in these files.

Code Overview

This module provides a high level model for CCD and CMOS staring array signal chain modelling. The work is based on a paper and Matlab code by Mikhail Konnik, available at:

See the documentation at http://nelisw.github.io/pyradi-docs/_build/html/index.html or pyradi/doc/rystare.rst for more detail.

Module functions

pyradi.rystare.photosensor(strh5, initialise=True)

This routine simulates the behaviour of a CCD/CMOS sensor, performing the conversion from irradiance to electrons, then volts, and then digital numbers.

The process from incident photons to the digital numbers appeared on the image is outlined. First of all, the radiometry is considered. Then, the process of conversion from photons to electrons is outlined. Following that, conversion from electrons to voltage is described. Finally, the ADC converts the voltage signal into digital numbers. The whole process is depicted on Figure below.

camerascheme-horiz.png

Many noise sources contribute to the resulting noise image that is produced by photosensors. Noise sources can be broadly classified as either fixed-pattern (time-invariant) or temporal (time-variant) noise. Fixed-pattern noise refers to any spatial pattern that does not change significantly from frame to frame. Temporal noise, on the other hand, changes from one frame to the next.

Note that in the sequence below we add signal and noise signals linearly together. For uncorrelated noise sources, the noise power values are added in quadrature, but that does not apply here, because we are adding instantaneous noise values (per pixel) so that these noise and signal values add linearly.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.set_photosensor_constants(strh5, initialise=True)

Defining the constants that are necessary for calculation of photon energy, dark current rate, etc.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.check_create_datasets(strh5, initialise=True)

Create the arrays to store the various image-sized variables.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.source_follower(strh5, initialise=True)

The amplification of the voltage from Sense Node by Source Follower.

Conventional sensor use a floating-diffusion sense node followed by a charge-to-voltage amplifier, such as a source follower.

_images/source_follower.png

Source follower is one of basic single-stage field effect transistor (FET) amplifier topologies that is typically used as a voltage buffer. In such a circuit, the gate terminal of the transistor serves as the input, the source is the output, and the drain is common to both input and output. At low frequencies, the source follower has voltage gain:

{A_{\text{v}}} = \frac{v_{\text{out}}}{v_{\text{in}}} = \frac{g_m R_{\text{S}}}{g_m R_{\text{S}} + 1} \approx 1 \qquad (g_m R_{\text{S}} \gg 1)

Source follower is a voltage follower, its gain is less than 1. Source followers are used to preserve the linear relationship between incident light, generated photoelectrons and the output voltage.

The V/V non-linearity affect shot noise (but does not affect FPN curve) and can cause some shot-noise probability density compression. The V/V non-linearity non-linearity is caused by non-linear response in ADC or source follower.

The V/V non-linearity can be simulated as a change in source follower gain A_{SF} as a linear function of signal:

A_{SF_{new}} = \alpha \cdot \frac{V_{REF} - S(V_{SF}) }{V_{REF} } + A_{SF},

where \alpha = A_{SF}\cdot\frac{\gamma_{nlr} -1}{ V_{FW} } and \gamma_{nlr} is a non-linearity ratio of A_{SF}. In the simulation we assume A_{SF} = 1 and \gamma_{nlr} = 1.05 i.e. 5% of non-linearity of A_{SF}. Then the voltage is multiplied on the new sense node gain A_{SF_{new}}:

I_{V} = I_{V}\cdot A_{SF_{new}}

After that, the voltage goes to ADC for quantisation to digital numbers.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.fixed_pattern_offset(strh5)

Add dark fixed pattens offset

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.cds(strh5)

Reducing the noise by Correlated Double Sampling, but right now the routine just adds the noise.

Correlated Double Sampling (CDS) is a technique for measuring photo voltage values that removes an undesired noise. The sensor’s output is measured twice. Correlated Double Sampling is used for compensation of Fixed pattern noise caused by dark current leakage, irregular pixel converters and the like. It appears on the same pixels at different times when images are taken. It can be suppressed with noise reduction and on-chip noise reduction technology. The main approach is CDS, having one light signal read by two circuits.

In CDS, a circuit measures the difference between the reset voltage and the signal voltage for each pixel, and assigns the resulting value of charge to the pixel. The additional step of measuring the output node reference voltage before each pixel charge is transferred makes it unnecessary to reset to the same level for each pixel.

First, only the noise is read. Next, it is read in combination with the light signal. When the noise component is subtracted from the combined signal, the fixed-pattern noise can be eliminated.

CDS is commonly used in image sensors to reduce FPN and reset noise. CDS only reduces offset FPN (gain FPN cannot be reduced using CDS). CDS in CCDs, PPS, and photogate APS, CDS reduces reset noise, in photodiode APS it increases it See Janesick’s book and especially El Gamal’s lectures.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.adc(strh5)

An analogue-to-digital converter (ADC) transforms a voltage signal into discrete codes.

An analogue-to-digital converter (ADC) transforms a voltage signal into discrete codes. An N-bit ADC has 2^N possible output codes with the difference between code being V_{ADC.REF}/2^N. The resolution of the ADC indicates the number of discrete values that can be produced over the range of analogue values and can be expressed as:

K_{ADC} = \frac{V_{ADC.REF} - V_\mathrm {min}}{N_{max}} where V_\mathrm{ADC.REF} is the maximum voltage that can be quantified, V_{min} is minimum quantifiable voltage, and N_{max} = 2^N is the number of voltage intervals. Therefore, the output of an ADC can be represented as:

ADC_{Code} = \textrm{round}\left( \frac{V_{input}-V_{min}}{K_{ADC}} \right)

The lower the reference voltage V_{ADC.REF}, the smaller the range of the voltages one can measure.

After the electron matrix has been converted to voltages, the sense node reset noise and offset FPN noise are added, the V/V gain non-linearity is applied (if desired), the ADC non-linearity is applied (if necessary). Finally the result is multiplied by ADC gain and rounded to produce the signal as a digital number:

I_{DN} =  \textrm{round} (A_{ADC}\cdot I_{total.V}),

where I_\textrm{total.V} = (V_{ADC.REF} - I_{V}) is the total voltage signal accumulated during one frame acquisition, V_{ADC.REF} is the maximum voltage that can be quantified by an ADC, and I_V is the total voltage signal accumulated by the end of the exposure (integration) time and conversion. Usually I_V = I_{SN.V} after the optional V/V non-linearity is applied. In this case, the conversion from voltages to digital signal is linear. The adcnonlinearity “non-linear ADC case is considered below”.

In terms of the ADC, the following non-linearity and noise should be considered for the simulations of the photosensors: Integral Linearity Error, Differential Linearity Error, quantisation error, and ADC offset.

The DLE indicates the deviation from the ideal 1 LSB (Least Significant Bit) step size of the analogue input signal corresponding to a code-to-code increment. Assume that the voltage that corresponds to a step of 1 LSB is V_{LSB}. In the ideal case, a change in the input voltage of V_{LSB} causes a change in the digital code of 1 LSB. If an input voltage that is more than V_{LSB} is required to change a digital code by 1 LSB, then the ADC has DLE error. In this case, the digital output remains constant when the input voltage changes from, for example, 2 V_{LSB} to 4 V_{LSB}, therefore corresponding the digital code can never appear at the output. That is, that code is missing.

_images/dle.png

In the illustration above, each input step should be precisely 1/8 of reference voltage. The first code transition from 000 to 001 is caused by an input change of 1 LSB as it should be. The second transition, from 001 to 010, has an input change that is 1.2 LSB, so is too large by 0.2 LSB. The input change for the third transition is exactly the right size. The digital output remains constant when the input voltage changes from 4 LSB to 5 LSB, therefore the code 101 can never appear at the output.

The ILE is the maximum deviation of the input/output characteristic from a straight line passed through its end points. For each voltage in the ADC input, there is a corresponding code at the ADC output. If an ADC transfer function is ideal, the steps are perfectly superimposed on a line. However, most real ADC’s exhibit deviation from the straight line, which can be expressed in percentage of the reference voltage or in LSBs. Therefore, ILE is a measure of the straightness of the transfer function and can be greater than the differential non-linearity. Taking the ILE into account is important because it cannot be calibrated out.

_images/ILE.png

For each voltage in the ADC input there is a corresponding word at the ADC output. If an ADC is ideal, the steps are perfectly superimposed on a line. But most of real ADC exhibit deviation from the straight line, which can be expressed in percentage of the reference voltage or in LSBs.

In our model, we simulate the Integral Linearity Error (ILE) of the ADC as a dependency of ADC gain A_{ADC.linear} on the signal value. Denote \gamma_{ADC.nonlin} as an ADC non-linearity ratio (e.g., \gamma_{ADC.nonlin}
= 1.04). The linear ADC gain can be calculated from Eq.~ref{eq:kadc} as A_{ADC} = 1/K_{ADC} and used as A_{ADC.linear}. The non-linearity coefficient \alpha_{ADC} is calculated as:

\alpha_{ADC} = \frac{1}{V_{ADC.REF}} \left( \frac{ \log(\gamma_{ADC.nonlin}
\cdot A_{ADC.linear} )}{\log(A_{ADC.linear})} - 1 \right)

where V_\mathrm{ADC.REF} is the maximum voltage that can be quantified by an ADC:

A_{ADC.nonlin} = A_{ADC.linear}^{1-\alpha_{ADC} I_{total.V}},

where A_{ADC.linear} is the linear ADC gain. The new non-linear ADC conversion gain A_{ADC.nonlin} is then used for the simulations.

Quantisation errors are caused by the rounding, since an ADC has a finite precision. The probability distribution of quantisation noise is generally assumed to be uniform. Hence we use the uniform distribution to model the rounding errors.

It is assumed that the quantisation error is uniformly distributed between -0.5 and +0.5 of the LSB and uncorrelated with the signal. Denote q_{ADC} the quantising step of the ADC. For the ideal DC, the quantisation noise is:

\sigma_{ADC} = \sqrt{ \frac{q_{ADC}^2 }{12}}.

If q_{ADC} = 1 then the quantisation noise is \sigma_{ADC} = 0.29 DN. The quantisation error has a uniform distribution. We do not assume any particular architecture of the ADC in our high-level sensor model. This routine performs analogue-to-digital convertation of volts to DN.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.charge_to_voltage(strh5)

The charge to voltage conversion occurs inside this routine

V/e nonlinearity is small for CCD detetors, but can be very high for some CMOS architectures (up to 200%) [from Janesick p87]

A new matrix strh5[‘rystare/signal/voltage’] is created and the raw voltage signal is stored.

After the charge is generated in the pixel by photo-effect, it is moved row-by-row to the sense amplifier that is separated from the pixels in case of CCD. The packets of charge are being shifted to the output sense node, where electrons are converted to voltage. The typical sense node region is presented on Figure below.

_images/CCD-sensenoderegion.png

The sense node is the final collecting point at the end of the horizontal register of the CCD sensor. The CCD pixels are made with MOS devices used as reverse biased capacitors. The charge is readout by a MOSFET based charge to voltage amplifier. The output voltage is inversely proportional to the sense node capacitor. Typical example is that the sense node capacitor of the order 50fF, which produces a gain of 3.2 \mu V/ e^-. It is also important to minimize the noise of the output amplifier, textbf{typically the largest noise source in the system}. Sense node converts charge to voltage with typical sensitivities 1\dots 4 \mu V/e^-.

The charge collected in each pixel of a sensor array is converted to voltage by sense capacitor and source-follower amplifier.

Reset noise is induced during such conversion. Prior to the measurement of each pixel’s charge, the CCD sense capacitor is reset to a reference level. Sense node converts charge to voltage with typical sensitivities 1\dots 4 \mu V/e^-. The charge collected in each pixel of a sensor array is converted to voltage by sense capacitor and source-follower amplifier. Reset noise is induced during such conversion. Prior to the measurement of each pixel’s charge, the CCD sense node capacitor is reset to a reference level.

Sense Node gain non-linearity, or V/e non-linearity

The V/e^- non-linearity affect both FPN and shot noise and can cause some shot-noise probability density compression. This type of non-linearity is due to sense node gain non-linearity. Then sense node sensitivity became non-linear (see Janesick’s book):

S_{SN} ( V_{SN}/e^- ) = \frac{S(V_{SN}) }{(k_1/q)  \ln( V_{REF}/[V_{REF} - S(V_{SN})] )}

The V/e^- non-linearity can be expressed as a non-linear dependency of signals in electron and a sense-node voltage:

S[e^-] = \frac{k1}{q} \ln \left[ \frac{V_{REF}}{ V_{REF} -  S(V_{SN}) } \right]

The V/e^- non-linearity affects photon shot noise and skews the distribution, however this is a minor effect. The V/e^- non-linearity can also be thought as a sense node capacitor non-linearity: when a small signal is measured, C_{SN} is fixed or changes negligible; on the other hand, C_{SN} changes significantly and that can affect the signal being measured.

For the simulation purpose, the V/e^- non-linearity can be expressed as:

V_{SN} = V_{REF} - S(V_{SN}) = V_{REF}\exp\left[ - \frac{\cdot S[e^-]\cdot q }{k1} \right]

where k1=10.909*10^{-15} and q is the charge of an electron. The nonlinear capacitance is given by C = k1/V

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.sense_node_reset_noise(strh5, seed=None)

This routine calculates the noise standard deviation for the sense node reset noise.

Sense node Reset noise (kTC noise)

Prior to the measurement of each pixel’s charge packet, the sense node capacitor is reset to a reference voltage level. Noise is generated at the sense node by an uncertainty in the reference voltage level due to thermal variations in the channel resistance of the MOSFET reset transistor. The reference level of the sense capacitor is therefore different from pixel to pixel.

Because reset noise can be significant (about 50 rms electrons), most high-performance photosensors incorporate a noise-reduction mechanism such as correlated double sampling (CDS).

kTC noise occurs in CMOS sensors, while for CCD sensors the sense node reset noise is removed~ (see Janesick’s book) by Correlated Double Sampling (CDS). Random fluctuations of charge on the sense node during the reset stage result in a corresponding photodiode reset voltage fluctuation. The sense node reset noise (in volt units) is given by:

\sigma_{RESET}=\sqrt{\frac{k_B T}{C_{SN}}}

By the relationship Q=CV it can be shown that the kTC noise can be expressed as electron count by

\sigma_{RESET}=\frac{\sqrt{k_B T C_{SN}}}{q} = \frac{k_B T}{q A_{SN}}

see also https://en.wikipedia.org/wiki/Johnson%E2%80%93Nyquist_noise

The simulation of the sense node reset noise may be performed as an addition of non-symmetric probability distribution to the reference voltage V_{REF}. However, the form of distribution depends on the sensor’s architecture and the reset technique. An Inverse-Gaussian distribution can be used for the simulation of kTC noise that corresponds to a hard reset technique in the CMOS sensor, and the Log-Normal distribution can be used for soft-reset technique. The sense node reset noise can be simulated for each (i,j)-th pixel for the soft-reset case as:

I_{SN.reset.V}=ln\mathcal{N}(0,\sigma_{RESET}^2)

then added to the matrix I_{REF.V} in Volts that corresponds to the reference voltage.

Note: For CCD, the sense node reset noise is entirely removed by CDS.

Note: In CMOS photosensors, it is difficult to remove the reset noise for the specific CMOS pixels architectures even after application of CDS. Specifically, the difficulties arise in ‘rolling shutter’ and ‘snap’ readout modes. The reset noise is increasing after CDS by a factor of \sqrt{2}. Elimination of reset noise in CMOS is quite challenging.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.dark_current_and_dark_noises(strh5)

This routine for adding dark current signals and noise, including dark FPN and dark shot noise.

This model is taken from Janesick’s ‘Photon Transfer’ book, page 168, which in turn is taken from Janesick’s ‘Scientific Charge-Coupled Devices’ book, page 622.

The dark signal is calculated for all pixels in the model. It is implemented using ones function in MATLAB as a matrix of the same size as the simulated photosensor. For each (i,j)-th pixel we have:

I_{dc.e^-} = t_I\cdot D_R,

where D_R is the average dark current (originally derived for silicon):

D_R = 2.55\cdot10^{15}P_A D_{FM} T^{1.5} \exp\left[-\frac{E_{gap}}{2\cdot k\cdot T}\right],

where: D_R is in units of [e-1/s], P_A is the pixel’s area [cm2]; D_{FM} is the dark current figure-of-merit in units of [nA/cm2] at 300K, varies significantly with detector material and sensor manufacturer, and used in this simulations as 0.5 nA/cm2 for silicon; E_{gap} is the bandgap energy of the semiconductor which also varies with temperature; k is Boltzman’s constant that is 8.617\cdot10^{-5} [eV/K].

The relationship between band gap energy and temperature can be described by Varshni’s empirical expression,

E_{gap}(T)=E_{gap}(0)-\frac{\alpha T^2}{T+\beta},

where E_{gap}(0), \alpha and \beta are material constants. The energy bandgap of semiconductors tends to decrease as the temperature is increased. This behaviour can be better understood if one considers that the inter-atomic spacing increases when the amplitude of the atomic vibrations increases due to the increased thermal energy. This effect is quantified by the linear expansion coefficient of a material.

For the Silicon: E_{gap}(0) = 1.1557 [eV], \alpha = 7.021*10^{-4} [eV/K], and
\beta = 1108 [K].
It appears that fill factor does not apply to dark noise (Janesick book p168 and Konnik’s code
does not show this).

According to Janesick’s Photon transfer book p169 the dark current FPN standard deviation is around 10% (CCD) and 40% (CMOS) of the dark current. Note that ‘dark’ FPN (DN) is much greater than ‘light’ FPN (PN) by approximately 10 to 40 times.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields, dark current in electrons
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.source_follower_noise(strh5, initialise=True)

The source follower noise routine, calculates noise in volts.

The pixel’s source follower noise limits the read noise, however in high-end CCD and CMOS cameras the source follower noise has been driven down to one electron rms. Pixel source follower MOSFET noise consists of three types of noise: - white noise; - flicker noise; - random telegraph noise (RTS). Each type of noise has its own physics that will be briefly sketched below.

Johnson noise (white noise)

Similarly to the reset noise in sense node, the source-follower amplifier MOSFET has a resistance that generates thermal noise whose value is governed by the Johnson white noise equation. It is therefore either referred to as Johnson noise or simply as white noise, since its magnitude is independent of frequency. If the effective resistance is considered to be the output impedance of the source-follower amplifier, the white noise, in volts, is determined by the following equation:

N_{white} (V_{SF}) = \sqrt{4kTBR_{SF}}

where k is Boltzmann’s constant (J/K), T is temperature [K], B refers to the noise power bandwidth [Hz], and R_{SF} is the output impedance of the source-follower amplifier.

Flicker noise

The flicker noise is commonly referred to as 1/f noise because of its approximate inverse dependence on frequency. For cameras in which pixels are read out at less than approximately 1 megahertz, and with a characteristic 1/f noise spectrum, the read noise floor is usually determined by 1/f noise. Note that the noise continues to decrease at this rate until it levels off, at a frequency referred to as the 1/f corner frequency. For the typical MOSFET amplifier, the white noise floor occurs at approximately 4.5 nV/Hz^{1/2}.

Prominent sources of 1/f noise in an image sensor are pink-coloured noise generated in the photo-diodes and the low-bandwidth analogue operation of MOS transistors due to imperfect contacts between two materials. Flicker noise is generally accepted to originate due to the existence of interface states in the image sensor silicon that turn on and off randomly according to different time constants. All systems exhibiting 1/f behaviour have a similar collection of randomly-switching states. In the MOSFET, the states are traps at the silicon-oxide interface, which arise because of disruptions in the silicon lattice at the surface. The level of 1/f noise in a CCD sensor depends on the pixel sampling rate and from certain crystallographic orientations of silicon wafer.

Random Telegraph Signal (RTS) noise

As the CCD and CMOS pixels are shrinking in dimensions, the low-frequency noise increases. In such devices, the low-frequency noise performance is dominated by Random Telegraph Signals (RTS) on top of the 1/f noise. The origin of such an RTS is attributed to the random trapping and de-trapping of mobile charge carriers in traps located in the oxide or at the interface. The RTS is observed in MOSFETs as a fluctuation in the drain current. A pure two-level RTS is represented in the frequency domain by a Lorentzian spectrum.

Mathematically the source follower’s noise power spectrum can be described as: S_{SF}(f) = W(f)^2 \cdot \left(1 + \frac{f_c}{f}\right)+S_{RTS}(f),

where W(f) is the thermal white noise [V/Hz^{1/2}, typically 15 nV/Hz^{1/2} ], flicker noise corner frequency f_c in [Hz] (flicker noise corner frequency is the frequency where power spectrum of white and flicker noise are equal), and the RTS power spectrum is given (see Janesick’s book):

S_{RTS}(f) = \frac{2\Delta I^2 \tau_{RTS}}{4+(2\pi f  \tau_{RTS})^2},

where \tau_{RTS} is the RTS characteristic time constant [sec] and \Delta I is the source follower current modulation induced by RTS [A].

The source follower noise can be approximated as:

\sigma_{SF} = \frac{\sqrt{\int\limits_{0}^{\infty} S_{SF}(f) H_{CDS}(f) df }}{A_{SN}A_{SF}(1-\exp^{-t_s/\tau_D})}

where: - \sigma_{SF} is the source follower noise [e- rms] - f is the electrical frequency [Hz] - t_s is the CDS sample-to-sampling time [sec] - \tau_D is the CDS dominant time constant (see Janesick’s Scientific CCDs book) usually set as \tau_D = 0.5t_s [sec].

The H_{CDS}(f) function is the CDS transfer function is (see Janesick’s book):

H_{CDS}(f) = \frac{1}{1+(2\pi f \tau_D)^2} \cdot [2-2\cos(2\pi f t_s)]

First term sets the CDS bandwidth for the white noise rejection before sampling takes place through B = 1/(4\tau_D), where B is defined as the noise equivalent bandwidth [Hz].

Note: In CCD photosensors, source follower noise is typically limited by the flicker noise.

Note: In CMOS photosensors, source follower noise is typically limited by the RTS noise. As a side note, such subtle kind of noises is visible only on high-end ADC like 16 bit and more.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.multiply_detector_area(strh5)

This routine multiplies detector area

The input to the model of the photosensor is assumed to be a matrix E_{q}\in R^{N\times M} that has been converted to electronrate irradiance, corresponding to electron rate [e/(m2.s)]. The electron rate irriance is converted to electron rate into the pixel by accounting for detector area:

\Phi_{q}  =  \textrm{round} \left(  E_{q} \cdot P_A    \right),

where P_A is the area of a pixel [m2].

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.multiply_integration_time(strh5)

This routine multiplies with integration time

The input to the model of the photosensor is assumed to be a matrix E_{q}\in R^{N\times M} that has been converted to electrons, corresponding to electron rate [e/s]. The electron rate is converted to electron count into the pixel by accounting for detector integration time:

\Phi_{q}  =  \textrm{round} \left(  E_{q} \cdot t_I  \right),

where t_{I} is integration (exposure) time.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.convert_to_electrons(strh5)

This routine converts photon rate irradiance to electron rate irradiance

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.shotnoise(sensor_signal_in)

This routine adds photon shot noise to the signal of the photosensor that is in photons.

The photon shot noise is due to the random arrival of photons and can be described by a Poisson process. Therefore, for each (i,j)-th element of the matrix \Phi_{q} that contains the number of collected photons, a photon shot noise is simulated as a Poisson process \mathcal{P} with mean \Lambda:

\Phi_{ph.shot}=\mathcal{P}(\Lambda), \,\,\,\,\mbox{ where   } \Lambda = \Phi_{q}.

We use the ryutils.poissonarray function that generates Poisson random numbers with mean \Lambda. That is, the number of collected photons in (i,j)-th pixel of the simulated photosensor in the matrix \Phi_{q} is used as the mean \Lambda for the generation of Poisson random numbers to simulate the photon shot noise. The input of the ryutils.poissonarray function will be the matrix \Phi_{q} that contains the number of collected photons. The output will be the matrix \Phi_{ph.shot} \rightarrow \Phi_{q}, i.e., the signal with added photon shot noise. The matrix \Phi_{ph.shot} is recalculated each time the simulations are started, which corresponds to the temporal nature of the photon shot noise.

Args:
sensor_signal_in (np.array[N,M]): photon irradiance in, in photons
Returns:
sensor_signal_out (np.array[N,M]): photon signal out, in photons
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.responsivity_FPN_light(strh5)

Multiploying the photon signal with the PRNU.

The Photo Response Non-Uniformity (PRNU) is the spatial variation in pixel conversion gain (from photons to electrons). When viewing a uniform scene the pixel signals will differ because of the PRNU, mainly due to variations in the individual pixel’s characteristics such as detector area and spectral response. These variations occur during the manufacture of the substrate and the detector device.

The PRNU is signal-dependent (proportional to the input signal) and is fixed-pattern (time-invariant). For visual (silicon) sensors the PRNU factor is typically 0.01\dots 0.05, but for HgCdTe sensors it can be as large as 0.02\dots 0.25. It varies from sensor to sensor, even within the same manufacturing batch.

The photo response non-uniformity (PRNU) is considered as a temporally-fixed light signal non-uniformity. The PRNU is modelled using a Gaussian distribution for each (i,j)-th pixel of the matrix I_{e^-}, as I_{PRNU.e^-}=I_{e^-}(1+\mathcal{N}(0,\sigma_{PRNU}^2)) where \sigma_{PRNU} is the PRNU factor value.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.responsivity_FPN_dark(strh5)

Add dark current noises that consist of Dark FPN and Dark shot noise.

Pixels in a hardware photosensor cannot be manufactured exactly the same from perfectly pure materials. There will always be variations in the photo detector area that are spatially uncorrelated, surface defects at the SiO_2/Si interface (see Sakaguchi paper on dark current reduction), and discrete randomly-distributed charge generation centres. These defects provide a mechanism for thermally-excited carriers to move between the valence and conduction bands. Consequently, the average dark signal is not uniform but has a spatially-random and fixed-pattern noise (FPN) structure. The dark current FPN can be expressed as follows:

\sigma_{d.FPN} = t_I D_R \cdot D_N,

where t_I is the integration time, D_R is the average dark current, and D_N is the dark current FPN factor that is typically 0.1\dots 0.4 for CCD and CMOS sensors.

There are also so called ‘outliers’ or ‘dark spikes’; that is, some pixels generate a dark signal values much higher than the mean value of the dark signal. The mechanism of such ‘dark spikes’ or ‘outliers’ can be described by the Poole-Frenkel effect (increase in emission rate from a defect in the presence of an electric field).

Simulation of dark current fixed pattern noise

The dark current Fixed Pattern Noise (FPN) is simulated using non-symmetric distributions to account for the ‘outliers’ or ‘hot pixels’. It is usually assumed that the dark current FPN can be described by Gaussian distribution. However, such an assumption provides a poor approximation of a complicated noise picture.

Studies show that a more adequate model of dark current FPN is to use non-symmetric probability distributions. The concept is to use two distributions to describe very ‘leaky’ pixels that exhibit higher noise level than others. The first distribution is used for the main body of the dark current FPN, with a uniform distribution superimposed to model ‘leaky’ pixels. For simulations at room-temperature (25^\circ C) authors use a logistic distribution, where a higher proportion of the population is distributed in the tails. For higher temperatures, inverse Gaussian and Log-Normal distributions have been proposed. The Log-Normal distribution works well for conventional 3T APS CMOS sensors with comparatively high dark current.

In our simulations we use the Log-Normal distribution for the simulation of dark current FPN in the case of short integration times, and superimposing other distributions for long integration times. The actual simulation code implements various models, including Log-Normal, Gaussian, and Wald distribution to emulate the dark current FPN noise for short- and long-term integration times.

The dark current FPN for each pixel of the matrix I_{dc.shot.e^-} is computed as:

I_{dc.FPN.e^-}  = I_{dc.shot.e^-}  + I_{dc.shot.e^-} \cdot ln\mathcal{N}(0,\sigma_{dc.FPN.e^-}^2)

where \sigma_{dc.FPN.e^-} = t_I D_R  D_N, D_R is the average dark current, and D_N is the dark current FPN factor. Since the dark current FPN does not change from one frame to the next, the matrix ln \mathcal{N} is calculated once and then can be re-used similar to the PRNU simulations.

The experimental results confirm that non-symmetric models, and in particular the Log-Normal distribution, adequately describe the dark current FPN in CMOS sensors, especially in the case of a long integration time (longer than 30-60 seconds). For long-exposure case, one needs to superimpose two (or more, depending on the sensor) probability distributions.

Args:
strh5 (hdf5 file): hdf5 file that defines all simulation parameters
Returns:
in strh5: (hdf5 file) updated data fields
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.FPN_models(sensor_signal_rows, sensor_signal_columns, noisetype, noisedistribution, spread, filter_params=None)

The routine contains various models on simulation of Fixed Pattern Noise.

There are many models for simulation of the FPN: some of the models are suitable for short-exposure time modelling (Gaussian), while other models are more suitable for log-exposure modelling of dark current FPN.

Gaussian model (Janesick-Gaussian)

Fixed-pattern noise (FPN) arises from changes in dark currents due to variations in pixel geometry during fabrication of the sensor. FPN increases exponentially with temperature and can be measured in dark conditions. Column FPN is caused by offset in the integrating amplifier, size variations in the integrating capacitor CF, channel charge injection from reset circuit. FPN components that are reduced by CDS. Dark current FPN can be expressed as:

\sigma_{D_{FPN}} = D\cdot D_N,

where D_N is the dark current FPN quality, which is typically between 10% and 40% for CCD and CMOS sensors (see Janesick’s book), and D = t_I D_R. There are other models of dark FPN, for instance as a autoregressive process.

El Gamal model of FPN with Autoregressive process

To capture the structure of FPN in a CMOS sensor we express F_{i,j} as the sum of a column FPN component Y_j and a pixel FPN component X_{i,j}. Thus, F_{i,j} = Y_j + X_{i,j}, where the Y_j’s and the X_{i,j}’s are zero mean random variables.

The first assumption is that the random processes Y_{j} and X_{i,j} are uncorrelated. This assumption is reasonable since the column and pixel FPN are caused by different device parameter variations. We further assume that the column (and pixel) FPN processes are isotropic.

The idea to use autoregressive processes to model FPN was proposed because their parameters can be easily and efficiently estimated from data. The simplest model, namely first order isotropic autoregressive processes is considered. This model can be extended to higher order models, however, the results suggest that additional model complexity may not be warranted.

The model assumes that the column FPN process Y_{j} is a first order isotropic autoregressive process of the form:

Y_j = a(Y_{j-1}+Y_{j+1}) + U_j

where the U_j s are zero mean, uncorrelated random variables with the same variance \sigma_U , and 0 \leq a \leq 1 is a parameter that characterises the dependency of Y_{j} on its two neighbours.

The model assumes that the pixel FPN process X_{i,j} is a two dimensional first order isotropic autoregressive process of the form:

X_{i,j} = b(X_{i-1,j} + X_{i+1,j} +  X_{i,j-1} + X_{i,j+1} ) + V_{i,j}

where the V_{i,j} s are zero mean uncorrelated random variables with the same variance \sigma_V , and 0 \leq b \leq 1 is a parameter that characterises the dependency of X_{i,j} on its four neighbours.

Args:
sensor_signal_rows(int): number of rows in the signal matrix
sensor_signal_columns(int): number of columns in the signal matrix
noisetype(string): type of noise to generate: [‘pixel’ or ‘column’]
noisedistribution(string): the probability distribution name [‘AR-ElGamal’, ‘Janesick-Gaussian’, ‘Wald’, ‘LogNormal’]
spread(float): spread around mean value (sigma/chi/lambda) for the probability distribution
filter_params([nd.array]): a vector of parameters for the probability filter
Returns:
noiseout (np.array[N,M]): generated noise of FPN.
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.nEcntLLightDF(tauAtmo, tauFilt, tauOpt, quantEff, rhoTarg, cosTarg, inttime, pfrac, detarea, fno, scenario, specBand, dfPhotRates)

Calculate the number of electrons in a detector given sensor parameters and photon radiance dataframe

All values in base SI units

Args:
tauAtmo (scalar or nd.array): atmosphere transmittance
tauFilt (scalar or nd.array): sensor filter transmittance
tauOpt (scalar or nd.array): sensor optics transmittance
quantEff (scalar or nd.array): sensor detector quantum efficiency
rhoTarg (scalar or nd.array): target diffuse reflectance
cosTarg (scalar): cos of illuminator angle wrt normal vector
inttime (scalar): integration time s
pfrac (scalar): fraction of clear optics
detarea (scalar): detector area m2
fno (scalar): f number
scenario (str): Scenario as key to rypflux.py dataframe
specBand (str): Spectral band as key to rypflux.py dataframe
dfPhotRadiance (pd.DataFrame): rypflux.py dataframe radiance in q/(s.m2.sr)
Returns:
n (float): number of electrons in charge well
Raises:
No exception is raised.

Author: CJ Willers

pyradi.rystare.nEcntLLightPhotL(tauAtmo, tauFilt, tauOpt, quantEff, rhoTarg, cosTarg, inttime, pfrac, detarea, fno, photRadiance)

Calculate the number of electrons in a detector given sensor parameters and photon radiance

All values in base SI units

Args:
tauAtmo (scalar or nd.array): atmosphere transmittance
tauFilt (scalar or nd.array): sensor filter transmittance
tauOpt (scalar or nd.array): sensor optics transmittance
quantEff (scalar or nd.array): sensor detector quantum efficiency
rhoTarg (scalar or nd.array): target diffuse reflectance
cosTarg (scalar): cos of illuminator angle wrt normal vector
inttime (scalar): integration time s
pfrac (scalar): fraction of clear optics
detarea (scalar): detector area m2
fno (scalar): f number
photRadiance (scalar): in-band photon radiance q/(s.m2.sr)
Returns:
n (float): number of electrons in charge well
Raises:
No exception is raised.

Author: CJ Willers

pyradi.rystare.nElecCntThermalScene(wl, tmptr, emis, tauAtmo, tauFilt, tauOpt, quantEff, inttime, pfrac, detarea, fno)

Calculate the number of electrons in a detector from a thermal source

All values in base SI units

Args:
wl (np.array): wavelength vector
tmptr (scalar): source temperature
emis (np.array of scalar): source emissivity
tauAtmo (scalar or nd.array): atmosphere transmittance
tauFilt (scalar or nd.array): sensor filter transmittance
tauOpt (scalar or nd.array): sensor optics transmittance
quantEff (scalar or nd.array): sensor detector quantum efficiency
inttime (scalar): integration time s
pfrac (scalar): fraction of clear optics
detarea (scalar): detector area m2
fno (scalar): f number
Returns:
n (float): number of electrons in charge well
Raises:
No exception is raised.

Author: CJ Willers

pyradi.rystare.nEcntThermalOptics(wl, tmptrOpt, tauFilt, tauOpt, quantEff, inttime, pfrac, detarea, fno)

Calculate the number of electrons in a detector from hot optics

All values in base SI units

Args:
wl (np.array): wavelength vector
tmptrOpt (scalar): optics temperature
tauFilt (scalar or nd.array): sensor filter transmittance
tauOpt (scalar or nd.array): sensor optics transmittance
quantEff (scalar or nd.array): sensor detector quantum efficiency
inttime (scalar): integration time s
pfrac (scalar): fraction of clear optics
detarea (scalar): detector area m2
fno (scalar): f number
Returns:
n (float): number of electrons in charge well
Raises:
No exception is raised.

Author: CJ Willers

pyradi.rystare.nElecCntReflSun(wl, tauSun, tauAtmo=1, tauFilt=1, tauOpt=1, quantEff=1, rhoTarg=1, cosTarg=1, inttime=1, pfrac=1, detarea=1, fno=0.8862269255, emissun=1.0, tmprt=6000.0)

Calculate the number of electrons in a detector or photon radiance for reflected sunlight

All values in base SI units.

By using the default values when calling the function the radiance at the source can be calculated.

Args:
wl (np.array (N,) or (N,1)): wavelength
tauSun (np.array (N,) or (N,1)): transmittance between the scene and sun
tauAtmo (np.array (N,) or (N,1)): transmittance between the scene and sensor
tauFilt (np.array (N,) or (N,1)): sensor filter transmittance
tauOpt (np.array (N,) or (N,1)): sensor optics transmittance
quantEff (np.array (N,) or (N,1)): detector quantum efficiency
rhoTarg (np.array (N,) or (N,1)): target diffuse surface reflectance
cosTarg (scalar): cosine between surface normal and sun/moon direction
inttime (scalar): detector integration time
pfrac (scalar): fraction of optics clear aperture
detarea (scalar): detector area
fno (scalar): optics fnumber
emissun (scalar): sun surface emissivity
tmprt (scalar): sun surface temperature
Returns:
n (scalar): number of electrons accumulated during integration time
Raises:
No exception is raised.
pyradi.rystare.darkcurrentnoise(inttime, detarea, temptr, Egap, DFM=5e-06)

Calculate the dark current noise given detector parameters

Args:
inttime (scalar): integration time in seconds
detarea (scalar): detector area in m2
temptr (scalar): temperature in K
Egap (scalar): bandgap in eV
DFM (scalar): in units of nA/m2
Returns:
n (scalar): dark current noise as number of electrons
Raises:
No exception is raised.
pyradi.rystare.kTCnoiseCsn(temptr, sensecapacity)
Args:
temptr (scalar): temperature in K
sensecapacity (): sense node capacitance F
Returns:
n (scalar): noise as number of electrons
Raises:
No exception is raised.
pyradi.rystare.kTCnoiseGv(temptr, gv)
Args:
temptr (scalar): temperature in K
gv (scalar): sense node gain V/e
Returns:
n (scalar): noise as number of electrons
Raises:
No exception is raised.
pyradi.rystare.define_metrics()

This simple routine defines various handy shorthand for cm and mm in the code.

The code defines a number of scaling factors to convert to metres and radians

Args:
None
Returns:
scaling factors.
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.limitzero(a, thr=0.6)

Performs an asymetric clipping to prevent negative values. The lower-end values are clumped up towards the lower positive values, while upper-end values are not affected.

This function is used to prevent negative random variables for wide sigma and low mean value, e.g., N(1,.5). If the random variables are passed through this function The resulting distribution is not normal any more, and has no known analytical form.

A threshold value of around 0.6 was found to work well for N(1,small) up to N(1,.5).

Before you use this function, first check the results using the code below in the main body of this file.

Args:
a (np.array): an array of floats,
Returns:
scaling factors.
Raises:
No exception is raised.

Author: CJ Willers

pyradi.rystare.run_example(doTest='Advanced', outfilename='Output', pathtoimage=None, doPlots=False, doHisto=False, doImages=False)

This code provides examples of use of the pyradi.rystare model for a CMOS/CCD photosensor.

Two models are provided ‘simple’ and ‘advanced’

doTest can be ‘Simple’ or ‘Advanced’

Args:
doTest (string): which example to run ‘Simple’, or ‘Advanced’
outfilename (string): filename for output files
pathtoimage (string): fully qualified path to where the image is located
doPlots (boolean): flag to control the creation of false colour image plots with colour bars
doHisto (boolean): flag to control the creation of image histogram plots
doImages (boolean): flag to control the creation of monochrome image plots
Returns:
hdffilename (string): output HDF filename
Raises:
No exception is raised.

Author: Mikhail V. Konnik, revised/ported by CJ Willers

Original source: http://arxiv.org/pdf/1412.4031.pdf

pyradi.rystare.get_summary_stats(hdffilename)

Return a string with all the summary input and results data.

Args:
hdffilename (string): filename for input HDF file
Returns:
Returns a string with summmary data.
Raises:
No exception is raised.

Author: CJ Willers

[hdf5asdataformatalt]https://github.com/NelisW/pyradi/blob/master/pyradi/hdf5-as-data-format.md