I'm currently working to simulate an image that mirrors the output of a CMOS camera we have in the lab. My images include Poisson (photon) noise as well as a gaussian (readout) noise.
My issue is dealing with the Poisson noise. I think I have an understanding here, but I've noticed a lack of consistency with the terminology (and sometimes implementation) so I was wondering if someone could clear this up for me, as much as I feel like a fool for asking.
For the rest of this, I will be using
image to reference my data array prior to the addition of noise. I will also use
poisson() to reference a function that takes an array of values as an input, and outputs an array of identical shape, having randomly generated new values using each elements input as the mean for a poisson distribution.
What I've Seen:
Poisson Noise as an additive to the image:
The application of this method makes the most sense to me as I generate images, as it's the general approach I've used for building images in the past, but I've only seen a couple implementations of this, and neither of their results make a lot of sense to me.
noisy_image = image + A*poisson(ones((len(image),len(image)))
noisy_image = image + poisson(image)
The first one has a totally arbitrary A and arbitrary input array. It is also not signal dependent, so I assume it's just flat out wrong.
The second one makes more sense since its output is actually signal dependent, but effectively doubles the values of the array. Surely making adjustments for noise should not double the values of my image?
Poisson Noise as an application of Poisson statistics to the expected count:
This kind of implementation was not the first to come to mind for me, though it is objectively just as, if not more, simple to do. That said, its results make a lot more sense to me as well.
noisy_image = poisson(image)
Each array value in this case is calculated with a Poisson distribution assuming the input value as the mean. This result is signal dependent, obeys the right distribution, and its output isn't total nonsense. That said, after applying a gain to the image, it actually comes out looking worse than the second additive listed above. That doesn't speak to its accuracy, but does make me want to confirm that my usage of it is correct.
This method also matches nicely with the (not especially helpful) statements claiming the size of the Poisson noise is sqrt(signal), as that is the standard deviation of this method.
So My Implementation is Then
If I'm ignoring the Gaussian noise for right this moment, my implementation has the form:
image = original_data_array
noisy_image = poisson(image)
final_output = (noisy_image*gain)+offset+background_light
This looks quite sensible to me, but poor wording in a few papers I've been looking at has me doubting myself. A confirmation of whether I'm looking at this sensibly or not would be of a great assistance.