exponential noise

 hi

   im new to image processing. i was trying to generate a matrix with random variables of exponential type(to model exponential noise). i used the cdf to generate them. 

 k=-1/a;
 R=k*log(1-rand(M,N));

my image is 256X256 sized, 0-255 grayvalued. problem is, no matter what values I give for 'a'(from .00001 to 10000), and compute R(which is my additive noise matrix), and then add it to my image I, I always get a full white image.i cant see even a speck of black?why is that? what values should i give?

Answer 1

[Walter Roberson]

Could you show the code you use?

log(1-rand(M,N)) would have rand() at least 2^(-53) and at most 1-2^(-53). 2^(-53) is more easily expresed as eps(1/2). Because of the 1 minus, the lowest rand() value translates to log(eps(1/2)) and the highest rand() value translates to log(1-eps(1/2)). Effectively the 1 minus step is wasted here.

log(eps(1/2)) is about -36 and log(1-eps(1/2)) is about -1E-16 . Your k is -1/a so the signs will be flipped, leading to one end being (1E-16)/a and the other end being 36/a . Which is larger will depend on whether a is less than 1 or greater. If a is greater than 72, then 36/a would be less than 1/2 (and (1E-16)/a much less still); when that is added to a uint8 value and the result converted back to a uint8 value, there would be no change. Unless, of course, you added it to the uint8 value an came up with a floating point result, in which case any value from 1 or higher would saturate as white. (hint hint)

If your a are much less than 1, even the 1E-16 would get magnified substantially. And if you have the double vs uint8 result as described above...

In short, you likely have a data type conversion problem.

[Walter Roberson]

Test out the class() of your input image, and the class() of your output image. What do you show as max(R(:)) ? How are you displaying R ?

[blackcoffee]

sorry.. i accidentally deleted what id posted.. so im reposting it

function R=createnoise(type, M, N, a, b) 
case 'exponential'
    if nargin<=3
        a=1;
    end
    if a<=0
        error('parameter a must be positive for exponential');
    end
    k=-1/a;
    R=k*log(1-rand(M,N));

..from Woods and GOnzalez.

both input and output image classes are char. I am not displaying R per se. Just adding it with my input image using imadd. And min value of R was 30. And Max 1.084e+04

[Walter Roberson]

Where is your switch() statement that goes with the 'case' ?

[Walter Roberson]

When you do the imadd(), is the image being added to uint8 data class? Is the result of imadd() still uint8 data class? What is the minimum and maximum values of the image after addition?

If the added image ended up being floating point, then as you have indicated that the minimum R is 30, every value after the addition would be at least 30, and so all of those floating point values would be greater than 1. The maximum value permitted for a floating point image is 1; everything beyond that is clamped to 1. The result would be an all-white image just like you are seeing, if you display the result using image() or imshow(). If you display the result using imagesc() you would see differences.

[blackcoffee]

switch lower(type)
        case 'uniform'
            R=a+(b-a)*rand([M,N]);
        case 'gaussian'
            R=a+b*randn(M,N);
        case 'lognormal'
            if nargin<=3
                a=1;b=.25;
            end
            R=a*exp(b*randn(M,N));
        case 'rayleigh'
            R=a+(-b*log(1-rand(M,N))).^.5;
        case 'exponential'
            if nargin<=3
                a=1;
            end
            if a<=0
                error('parameter a must be positive for exponential');
            end
            k=-1/a;
            R=k*log(1-rand(M,N));
        case 'salt & pepper'
            if nargin<=3
                a=.05;b=.05;
            end
          R(1:M,1:N) = .5;
            X(1:M,1:N)=rand(M,N);
            d=find (X<=a);
            R(d)=0;
            c=a+b;
            d=find (X>a & X<=c);
            R(d)=1;
        case 'erlang'
            if nargin<=3
                a=2;b=5;
            end
            if (b~=round(b)| b<=0)
                error('b must be a positive integer for erlang');
            end
            k=-1/a;
            R=zeros(M,N);
            for j=1:b
            R=R+k*log(1-rand(M,N));
            end
        otherwise
            error('unknown distribution type');
    end

that is my case().. m and n are sizes of the input image and a and b are probability or other parameters. the result of my imadd() is a float.. and when i use the imagesc function, for some unknown reason, the image doesn't show.. all i can see is big speckles of gray and black :( many thanks for your help so far.