2

How do I add Poisson noise with a specific mean of Poisson distribution in Matlab? The code below gives me a completely white image. 

img = im2double(phantom(512));
figure; imshow(img)
mean = 1e9;
R = poissrnd(mean, [512, 512]);
J = img + R;
figure; imshow(J);
Cris Luengo
  • 55,762
  • 10
  • 62
  • 120
ddxxdds
  • 373
  • 5
  • 13

2 Answers2

4

Typically, when adding Poisson noise to an image, you want to use the pixel's value as the mean (after scaling appropriately). This more closely would match Poisson noise in image acquisition.

For example, say your image, which is in the range 0-1 after im2double, was acquired in a photon-limited microscope. And let's say that the scaling is such that 1e9 photons on one pixel are represented as a pixel value of 1. To simulate the photon noise here, this one pixel with a (noiseless) value of 1 will have a (noisy) value of poissrnd(1e9)/1e9. Another pixel with a value of 0.5 will have a value of poissrnd(1e9*0.5)/1e9.

In the Image Processing Toolbox there is a function imnoise that does exactly this when called as follows:

img = imnoise(img,'poisson');

This function uses a scaling of 1e12. You wanted a scaling of 1e9, so you'll have to divide your image by 1e3:

img = imnoise(img/1e3,'poisson')*1e3;
Cris Luengo
  • 55,762
  • 10
  • 62
  • 120
2

The way you are doing it is good. Your mean is 1000000000 therefore the numbers are very very big.

However to show an image with auto-scaled colormap, do:

img = im2double(phantom(512));
figure; imshow(img,[])
mean = 1e9;
R = poissrnd(mean, [512, 512]);
J = img + R;
figure; imshow(J,[]);

Note the square brakets.

Ander Biguri
  • 35,140
  • 11
  • 74
  • 120
  • Note: I think you are just plainly doing whatever you mean to do wrong. Poisson distribution is discrete, i.e. it will give you integers, 0,1,2,.... However your image is [0-1], so just with `mean=1`, the noisy image will have a scale of ~x10 the initial value. You may want to add your noise as: `J=img+R*0.1*img;` or something like that. – Ander Biguri Aug 31 '18 at 16:30