PIC Logarithms

Here's an 8-bit logarithm algorithm that executes between 15 and 103 cycles. This version is based on a table look-up plus first order interpolation. See the logarithm theory page for more details.

        list    p=16C64,t=ON,c=132,n=80,st=off
        radix   dec


        include "P16C64.INC"

  cblock  0x20
        x,q,gi,gf,j,d

        z,xt

        interp,l1,l2x
        xl,xh,xstep
  endc



        ORG     0               ;Reset Vector

        GOTO    Main

        ORG     4               ;Interrupt Vector


Main

        BCF     STATUS,RP0      ;Point to BANK 0

        CLRF    xl
        CLRF    xh
        MOVLW   1
        MOVWF   xstep

xxx
        MOVF    xl,W
        MOVWF   x

        NOP
        CALL    log2_interp
        NOP

        MOVF    xstep,W
        ADDWF   xl,F
        SKPNC
          INCF  xh,F

        GOTO    xxx
;-------------------------------------------
;log2_interp
;
; The purpose of this routine is to find the base 2 logarithm of an
;8-bit integer:
;    gi.gf = lg(x)
; where,
;  x is the 8 bit input
;  lg() is the notation for base-2 logarithms
;  gi is interger portion of the result ( 0 <= gi <= 8) 
;  gf is the fractional portion of the result
;
;  The input variable, x, is first scaled to a value between 1 and 2. Then
;the fractional portion of the scaled x is multiplied by 256. The upper
;nibble is used as an index into the 16-entry logarithm table, while the
;lower nibble interpolates between successive entries in the table. The
;log table contains:
;  log_table[i] = 256 * (lg(1 + i/16))
;In other words, it contains the 16 equally spaced steps between 1 and 2.
;The multiplication by 256 scales the fraction to an 8 bit integer.
;
;psuedo code:
;
;  char x,i,j,gi,gf;
;  float xt;
;
;  xt = x/256;    // Actually this is an "implicit" division in assembly
;  gi = 0;
;  while(xt<1) {
;    xt = xt*2;      // Keep on shifting until xt is greater than 1.
;    gi = gi + 1;    // Count the number of shifts.
;  }
;
;  i = (xt-1)*256;   // Convert fractional portion to an integer.
;  j = i >> 4;
;  i = i & 0xf
;  gf = log_table[j] + i * (log_table[j+1] - log_table[j])
;
; inputs: x
; outputs: gi,gf
; RAM: l1,interp
;
; Execution time: 108 - 6*gi cycles if gf != 0
;                  60 - 6*gi cycles if gf == 0
; Worst case (when x=3) 103 cycles
; Best case (when x=128) 15 cycles

log2_interp

        MOVLW   HIGH(log_table_interp)
        MOVWF   PCLATH

        MOVLW   8               ;assume x=2^8
        MOVWF   gi
        CLRF    gf

        CLRC
lgi1    RLF     x,F             ;Search for the MSB of x
        SKPNC
          goto  lgi2
        DECFSZ  gi,F
          goto  lgi1

    ;If we get here, then there is an error (i.e. x==0)
        DECF    gf,F

lgi2    DECF    gi,F
        MOVF    x,F             ;If x is a perfect power of 2, then we're done
        SKPNZ
          RETURN

    ;Get the 3-bit interpolation factor that is in the lower nibble of x.
    ;The reason for three bits is because the 8 bit input was shifted left
    ;at least once leaving only 7 bits. The upper 4 are used for the
    ;index, while the lower three are used for the interpolation factor.

        RRF     x,W
        ANDLW   0x7
        IORLW   0x80            ;Set the MSB. Used below in the multiply loop
                                ;as the "loop counter".
        MOVWF   interp

        SWAPF   x,W
        ANDLW   0xf
        ADDLW   1
        MOVWF   gf                      ;Temporarily store the index
        CALL    log_table_interp        ;Get l2=log2( (x>>4) + 1)
        MOVWF   l1                      ;Store temporarily in l1

        DECF    gf,W                    ;
        CALL    log_table_interp        ;Get l1=log2( (x>>4) )

        SUBWF   l1,W                    ;W=l2-l1, This is always positive. 
        SUBWF   l1,F                    ;l1 = l1 - (l1-W) = W

    ;Initialize the product to zero
        CLRF    gf

    ;multiply interp and l2 - l1 and divide by 8. This is actually a 3 by 8
    ;bit multiplication. The division by 8 is implemented with a shift right
    ;one position for each of the three passes through the loop. Three passes 
    ;are determined by monitoring the migration of the MSB of interp. It is
    ;initialized to 1 above and is shifted right one position each pass
    ;through the loop. On the third pass, it will be in bit position 4.

lgi3    RRF     interp,F        ;If the LSB is set
        SKPNC
          ADDWF gf,F            ;Then add (l2-l1) to the product
        RRF     gf,F            ;Divide the product by two
        BTFSS   interp,4        ;If the MSB has not made it here yet
          goto  lgi3            ;then we need to loop again

        MOVF    l1,W
        ADDWF   gf,F
        RETURN

log_table_interp
        ADDWF   PCL,F

        RETLW   0x00    ;0x100*log2(0x100/0x100)
        RETLW   0x16    ;0x100*log2(0x110/0x100)
        RETLW   0x2b    ;0x100*log2(0x120/0x100)
        RETLW   0x3f    ;0x100*log2(0x130/0x100)
        RETLW   0x52    ;0x100*log2(0x140/0x100)
        RETLW   0x64    ;0x100*log2(0x150/0x100)
        RETLW   0x76    ;0x100*log2(0x160/0x100)
        RETLW   0x86    ;0x100*log2(0x170/0x100)
        RETLW   0x96    ;0x100*log2(0x180/0x100)
        RETLW   0xa5    ;0x100*log2(0x190/0x100)
        RETLW   0xb3    ;0x100*log2(0x1a0/0x100)
        RETLW   0xc1    ;0x100*log2(0x1b0/0x100)
        RETLW   0xcf    ;0x100*log2(0x1c0/0x100)
        RETLW   0xdc    ;0x100*log2(0x1d0/0x100)
        RETLW   0xe8    ;0x100*log2(0x1e0/0x100)
        RETLW   0xf4    ;0x100*log2(0x1f0/0x100)
        RETLW   0x00    ;0x100*log2(0x200/0x100)

        END

And here's more software.

BACK HOME



This page is maintained by Scott Dattalo. You can reach me at: scott@dattalo.com.
Last modified on 13DEC99.