Matlab digits precision ruins the calculation. How to avoid it?

Question

I work on an image, and i normalized it by its maximum value (3165 in this case).

After some work, I wanted to do some operations with the original values of the pixels, so I multiplied back with 3165.

The problem is that because of the digits precision, not all the pixels get back their integer value. For example consider a pixel with value of 30:

a=30/3165

a =

    0.0095

>> b=a*3165

b =

   30.0000

>> b-30

ans =

   3.5527e-15

The error itself is very small, but i need to use it as an index for other vector, and it's not work, obviously.

Technically, i can use round after multiplying. It's a good solution, but ugly, and not feels like the appropriate way and good practice to solve this issue.

Someone has any other (more professional) ideas? Or this is the way to do it?

Thanks.

Answer 1

The techniques mentioned below use eps and fix. These might be less "uglier" than rounding for your case. Let's look at these step by step with some codes and comments.

%%// Original indices
org_indices = [15 30 45 60]

%%// Scaled indices
scaled_indices = org_indices/3165

%%// Retrieved indices that are yet to be converted to exact integers
retrieved_indices = scaled_indices*3165

%%// Verify precision issue
check_ind1 = retrieved_indices(1) - 15

check_ind1 outputs 1.7764e-15 and is not exactly zero.

We can use eps that according to the official MATLAB documentation says - If d = eps(X), then d is the positive distance from abs(X) to the next larger in magnitude floating point number of the same precision as X.

Thus, we can calculate such distances/offsets as follows -

eps_val = eps(retrieved_indices)

eps_val could be thought of as window vector that determines within which if an element lies, must correspond to the integer that also lies within it. To make it easy for us, we can use the maximum of these window values for our decision making and to make it very safe, we can increase this width 10 times and let's call it max_window_width -

max_window_width = 10*max(eps_val)

This is not as ugly as rounding because rounding uses a window length of 1 because of its range [-0.5 0.5], whereas this one uses a window length of 7.1054e-14.

Finally, we can get the original indices back by adding this maximum width and then truncating the decimal portion using fix -

retrieved_indices2 = fix(retrieved_indices + max_window_width)

One can check back the values for exactness -

check_ind2 = retrieved_indices2(1) - 15

check_ind2 outputs 0.