EmbeddedRelated.com
Forums

Simple/Fast Algorithm for binary logarithm in C (on ARM7)

Started by Tilmann Reh May 11, 2007
"Spoon" <devnull@localhost.com> wrote in message 
news:4644647e$0$24673$426a34cc@news.free.fr...
> Peter Dickerson wrote: > >> I had a need for log of an integer with the result fixed point. For fixed >> point input there is just a constant difference. I didn't hugely accurate >> results, but the results are better than I needed. Input was signed >> 32-bit integer and result was ln(x) as 32-bit fixed point with 24-bit >> fraction. log10(x) was then obtained by scaling. There is no particular >> reason why the input can't be unsigned. > > Since the logarithm function is undefined for negative numbers, > why would an implementation allow anything but unsigned input?
Log is also undefined for 0 either. My problem was that my input could be negative (difference of a light signal from a dark signal) and so was easy to handle that way i.e propagate the error though the log. Peter
On May 11, 7:14 pm, Tilmann Reh <tilmann...@despammed.com> wrote:
> Hello all, > > in a current ARM7 project I need to do some calculations which include a > logarithm. All other math is done with integers (i.e. fixed point > arithmetic), and I would like to calculate the log fast and with little > code - so I don't prefer using FP conversion and the standard math > libraries. > > Can anyone point to fixed-point logarithmic routines I could use (resp. > tailor to my needs)? > > Thanks, > Tilmann > > --http://www.autometer.de- Elektronik nach Ma=DF.
Tilmann Wrote:
>=2Ein a current ARM7 project I need to do some calculations which include =
a logarithm. Sounds comparison operation to me as you need to compare the values with your logarithm table. IMHO hash table would yield the best performance with index and key attributes. ali
"Tilmann Reh" <tilmannreh@despammed.com> wrote in message 
news:464460b3$0$6401$9b4e6d93@newsspool2.arcor-online.net...
> Peter Dickerson schrieb: > >> I had a need for log of an integer with the result fixed point. For fixed >> point input there is just a constant difference. I didn't hugely accurate >> results, but the results are better than I needed. Input was signed >> 32-bit >> integer and result was ln(x) as 32-bit fixed point with 24-bit fraction. >> log10(x) was then obtained by scaling. There is no particular reason why >> the >> input can't be unsigned. >> >> However, the method requires a 256-entry 32-bit int table, a 256-entry >> byte >> table and one integer divide :-( >> >> I can't give out the source code because of copyright constraints but I'm >> happy to explain the algorithm. > > So, please start the explanation. :-)
First we shift the input left until a one is shifted out while subtracting the fixed point representation of ln(2) from the answer - initialised to 32*ln(2). The typical normization step. Obviously this stage can be speeded up with binary search approaches. So now we have the value 1+x+y in the form of a 32-bit fixed point number x+y with the binary point at the left. x, here, is the top eight bits and y the next 24. So 0 <= x < 1 is teps of 1/256 and 0 <= y < 1/256. We want ln(1+x+y) = ln((1+x)*(1+y/(1+x)) = ln(1+x) + ln(1+y/(1+x)). (a) The first term of (a) has only 256 values, and so can be looked up from a table. The second term is small, which we'll now estimate. Define z = y/(1+x), then ln(1+z) is approximately z - z*z/2 +... not that the z cubed and higher terms are smaller than 1 lsb and to a good approximation we can correct ln(1+z) from the linear with a small table because no more than eight bits are significant. OK, to make this a bit more explicit, lets use upper case letters for the binary (whole) numbers representing the various fractions. X = x*(2**8); Y = y*(2**32). Then the original normalized number is 1 + X/(2**8) + Y/(2**32). Look up ln(1+X/(2**8)) in a table indexed by X using a 24-bit fraction format - call this L. Set Z = Y/((2**8)+X) = z * (2**24), so its already in the 24-bit fraction format and less than 2**16. Add Z to L. Subtract a correction corresponding to the z*z/2 term by using (Z>>8) as an index into a byte-wide table. In the 24-fraction format the correction is roughly Z*Z/2/(2**24) or (Z/(2**8))*(Z/(2**8))/(2**9). So max correction is <128. So the longest step is the divide, which is actually 24-bit/9-bit. This could be speeded up on ARM7 by using the 32x32->64 multiply, looking up the multipliers corresponding to 1/(1+x). Hope this helps. I can explain further if necessary. Peter

Tilmann Reh wrote:

> in a current ARM7 project I need to do some calculations which include a > logarithm. All other math is done with integers (i.e. fixed point > arithmetic), and I would like to calculate the log fast and with little > code - so I don't prefer using FP conversion and the standard math > libraries. > > Can anyone point to fixed-point logarithmic routines I could use (resp. > tailor to my needs)?
It all depends on how accurate the log function should be. If you are fine with 10% accuracy, just count the insignificant bits and then do a simplest linear approximation. Vladimir Vassilevsky DSP and Mixed Signal Design Consultant http://www.abvolt.com
"Vladimir Vassilevsky" <antispam_bogus@hotmail.com> wrote in message 
news:us%0i.2641$zj3.1035@newssvr23.news.prodigy.net...
> > > Tilmann Reh wrote: > >> in a current ARM7 project I need to do some calculations which include a >> logarithm. All other math is done with integers (i.e. fixed point >> arithmetic), and I would like to calculate the log fast and with little >> code - so I don't prefer using FP conversion and the standard math >> libraries. >> >> Can anyone point to fixed-point logarithmic routines I could use (resp. >> tailor to my needs)? > > > It all depends on how accurate the log function should be. If you are fine > with 10% accuracy, just count the insignificant bits and then do a > simplest linear approximation.
The OP didn't give an input range or accuracy requirement so its hard to guess what is best. The thing we do know is that its for an ARM7, which has a reasonable multiplier, so better than linear shouldn't be a issue. Peter
On Fri, 11 May 2007 11:14:14 +0200, Tilmann Reh
<tilmannreh@despammed.com> wrote:

>in a current ARM7 project I need to do some calculations which include a >logarithm. All other math is done with integers (i.e. fixed point >arithmetic), and I would like to calculate the log fast and with little >code - so I don't prefer using FP conversion and the standard math >libraries. > >Can anyone point to fixed-point logarithmic routines I could use (resp. >tailor to my needs)?
You've gotten some responses already, but here's another you might consider looking at. It's one of the chapters of a book that Analog used to give away (don't know if they still do) to customers. Great book. Anyway, this chapter covers the logarithm for the ADSP-21xx DSP processor, but it includes a good description about the whys and wherefores as well as a coded example (which requires a little bit of learning of the isntruction set if you want to understand each and every step of the code -- but that isn't too hard.) Here is the chapter: http://www.analog.com/UploadedFiles/Associated_Docs/727652249Chapter_4.pdf Also, if you do actually decide to dig more deeply into the algorithm described and then coded up in the above chapter, you could read up on the instructions with this ZIP: http://www.analog.com/UploadedFiles/Associated_Docs/160158757adsp_21xx_um_rev3_0.zip (I couldn't find just the chapters online, quickly. So that is the whole thing, I imagine.) Jon
Rafael Deliano wrote:
>> Can anyone point to fixed-point logarithmic routines I could use >> (resp. tailor to my needs)? > > Good book with an idea on that : > Crenshaw "Math Toolkit for Real Time Programming" CMP Books 2000 >
Jack Crenshaw explains the whole thing (sprinkled throughout the book, but the area that especially applies here is on page 324, equation 10.53). I've used it to write fast approximation functions for several projects of the years. Log in a nutshell: Use log2() and scale the output as needed. The highest set bit tells you the left-of-radix value. The remaining bits represent the input in a range from 1.0 to 1.9999... Choose a table size (2^N) and the number of terms in the interpolation polynomial (3 usually works very well). Generate an A[] table of the log2(x) where x goes from 1.0 to 1+((2^N-1)/(2^N)). Make a table of "forward differences" (call it F) where F[i]=A[i+1]-A[i]. Make another difference table, G[i]=F[i+1]-F[i]. Now you can make B[] and C[] tables for the run-time code. C[i]=G[i]/2 B[i]=F[i]-G[i] run time code is: normalize x and remove MSB int i = x>>alphaBits // get table index int alpha = ((1<<alphaBits)-1) & x; // get delta between 1.0 and 1.999 retval = (C[i] * alpha) >> shift1 retval += B[i ] retval *= alpha retval >>= shift2 retval += A[i] The shift1 and shift2 are a function of your table size (N) and any extra scaling you might have done to improve the resolution of B[] and C[]. Crenshaw explains it all with Taylor series and shows how easy it is to extend to D[] or higher order tables. It comes down to a trade-off between accuracy, speed, and table space. I usually start by analyzing the error vs. table size in a spreadsheet, then use a C program so that I can test it over all possible input values. The C code might get hand converted to asm if I need the speed. HTH, Bob
Tilmann Reh wrote:
> > in a current ARM7 project I need to do some calculations which > include a logarithm. All other math is done with integers (i.e. > fixed point arithmetic), and I would like to calculate the log > fast and with little code - so I don't prefer using FP conversion > and the standard math libraries. > > Can anyone point to fixed-point logarithmic routines I could use > (resp. tailor to my needs)?
If you are using integers the logs are also integers. The bit position of the most significant bit gives the binary log. If you wish you can also round it up if the next two bits are 11. This is probably about as good as you can do in integers. You can extract this with (untested) while (b = c & (c - 1)) do c = b; /* now c is the log2 of the original c value */ It may be worthwhile to test for a zero value before starting this. -- <http://www.cs.auckland.ac.nz/~pgut001/pubs/vista_cost.txt> <http://www.securityfocus.com/columnists/423> <http://www.aaxnet.com/editor/edit043.html> <http://kadaitcha.cx/vista/dogsbreakfast/index.html> cbfalconer at maineline dot net -- Posted via a free Usenet account from http://www.teranews.com
Peter Dickerson schrieb:

> [Good explanation of log calculation]
Thank you very much, Peter, for that detailed explanation. Meanwhile, I also had a closer look at a combination of table look up for say the higher 8 bits of the mantissa, combined with simple linear interpolation for the rest. Doing that all in binary (with the binary log as result), this seems to be very efficient since no division is needed (can be substituted with shifting, which is an advantage especially with ARM7). The resulting precision also is very good (i.e. good enough for this application) with this method. If not, I could simply increase the table up to some reasonable amount since flash memory is not an issue here... However, I will archive your description for further reference. Thanks again for your help. Tilmann -- http://www.autometer.de - Elektronik nach Ma&#4294967295;.
Peter Dickerson schrieb:

>> It all depends on how accurate the log function should be. If you are fine >> with 10% accuracy, just count the insignificant bits and then do a >> simplest linear approximation. > > The OP didn't give an input range or accuracy requirement so its hard to > guess what is best. The thing we do know is that its for an ARM7, which has > a reasonable multiplier, so better than linear shouldn't be a issue.
10% will definitely not be enough... The input data is in fixed point format, so the range can be choosen rather freely - however the differences in the input data already are small (say about 1/1000), and the flatness of the log curve makes them even smaller in the log result. So the mantissa surely has to be taken into account... Tilmann -- http://www.autometer.de - Elektronik nach Ma&#4294967295;.