|
|
You can scale any positive real number by a power of 4 so that it lies in the interval [0.25,1).
Assume you're working with 32 bit integers that represent "0.32" fixed point numbers. All you need to do is shift left by two (=multiply by 4) until the value is >=0x40000000. That is the "pre-scaling" part. You count how many such shifts you needed (this will obviously always be less than 16, in the case of 32 bit numbers), call that number N. Then you calculate the square root. To "post-scale" (that is, undo the prescaling), you now have to divide by sqrt(4)^N=2^N, which is just a shift right by N bits.
In practice you can do a bit better than shifting up to 15 times to determine N (and normalize the input), because you don't have to do a "linear search" for the right scale factor; you can determine N with a simple divide-and-conquer-technique, assuming that x != 0:
int N = 0, t = x;
if(t & 0xffff0000) t >>= 16; else N += 8;
if(t & 0x0000ff00) t >>= 8; else N += 4;
if(t & 0x000000f0) t >>= 4; else N += 2;
if((t & 0x0000000c) == 0) N++;
|
Understanding that algorithm is left as an excercise to the reader :) - but seriously, this is a lot harder to write down than to understand; the general idea is that we want to count pairs of zero bits starting from the most significant bit of our number. If the first 16 bits are zero, we have 8 such pairs in the most significant 16 bits, plus whatever leading "zero packs" we have in the least significant 16 bits. If they're not zero, then the number of leading "zero packs" is the same as in the number (x >> 16), so in both cases we need to determine the number of leading zeroes of a 16-bit number next, which can be reduced exactly the same way.
You then need to shift x left by 2N bits before the sqrt computation and right by N bits afterwards. With fixed point it's a bit trickier; normally the easiest solution is just treating each fixed point number as if it was 0.32 then compensating for the actual number of fractional bits together with the "post-adjustment" shift.
|