LinuxQuestions.org - How to write a fixed point function?

Perhaps an illustrated example will help you understand.

Consider the linear blending function

Code:

float blend(const float x, const float y0, const float y1)

{

    return (1.0f - x) * y0 + x * y1;

}

The parameter x is between 0 and 1, inclusive. The function simply returns the linearly interpolated value y(x) where y(0)=y0 and y(1)=y1.

To use fixed-point math, you scale the floating-point value by a fixed integer S . So, instead of 0 to 1, integer x would get integer values from 0 to S .

When multiplying two fixed-point values, you must divide the result by S, because

Code:

S * x * y = (S * x) * (S * y) / S

With integer computation, you wish to do the division last, because intermediate results are also always integers, i.e.

Code:

(S * x) * (S * y) / S = trunc( ( trunc(S * x) * trunc(S * y) ) / S ) with integers

For example,

Code:

( ( x * y ) / 256 ) * 256 * 256 ) = trunc(x * y / 256) * 256 * 256 != 256 * x * y with integers

For similar reasons, you need to multiply a value by S before dividing it by a fixed-point value. Addition and subtraction need no special handling.

Finally, this code may not work with optimizations disabled:

Code:

avgRTT = ((255/256)*avgRtt) + ((1/256)*rtt);

because the divisions evaluate to zero; it is equivalent to

Code:

avgRTT = (trunc(255/256) * avgRtt) + (trunc(1/256) * rtt);

and evaluates to zero. However, if you are using any sort of optimization, the compiler masks the issue by reordering the operations.

Without explicit casts, C compilers may always reorder the order of the operations, as long as no precision is lost. They will never evaluate x*y/256 as trunc(x/256)*y . When optimizations are disabled, the C compiler may not reorder the operations at all. Explicit casts (in C99) require the compiler to evaluate the expression so that the cast part is limited to the precision of the cast type; i.e. (int)(x/256)*y will always be equivalent to trunc(x/256)*y in C99, regardless of optimization flags.

When the multiplication and division are done in the correct order,

Code:

avgRTT = (255 * avgRtt + rtt) / 256;

the code always works. That is, it works if and only if 255*avgRtt and 255*avgRtt+rtt fit into the integer type used; otherwise the result is garbage. It is evaluated effectively as

Code:

avgRtt = trunc(trunc(trunc(255 * avgRtt) + rtt) / 256);

In other words, when you work with integers, you must remember the implicit truncation to integers at every single step. True, it is a bit complicated, because without explicit casts the C compiler is allowed to optimize the order of the operations (as long as mathematical equivalence remains), so you cannot always tell where the implicit truncation would go.. but the reordered expressions never lose precision compared to the original: they may have more, but never less precision.

Returning back to the initial blending function: Assume x varies between 0 and 255, with scale factor S = 256, i.e. 256 = 1.0, and that we assume both y0 and y1 are integers that can be multiplied by 256 without overflowing. Then,

Code:

int blend(const int x, const int y0, const int y1)

{

    return ((256 - x) * y0 + x * y1) / 256;

}

is the equivalent fixed-point function.

Finally, we can optimize away one multiplication in both floating-point and integer cases. The end result is actually more numerically stable (for floats), too:

Code:

float blendf(const float x, const float y0, const float y1) { return y0 + x * (y1 - y0); }

int blendi(const int x, const int y0, const int y1) { return y0 + (x * (y1 - y0)) / 256; }

This function can even be written for any integer scale factor S, simply by supplying the scaled value of 1.0 to the function too:

Code:

int blend(const int x, const int one, const int y0, const int [I]y1[I])

{

    return (y0 + x * (y1 - y0)) / one;

}

Note, however, that if the intermediate results overflow (do not fit in an int), the result will be garbage. That is always an issue with fixed-point math.

Also note that all the variables above are either floats or (scaled) integers; I am not using any special fixed-point types at all. I am just using fixed-point arithmetic.

Did this help?