chiark - git - mdw - disorder/blob - lib/bits.c

   1 /*
   2  * This file is part of DisOrder
   3  * Copyright (C) 2008 Richard Kettlewell
   4  *
   5  * This program is free software; you can redistribute it and/or modify
   6  * it under the terms of the GNU General Public License as published by
   7  * the Free Software Foundation; either version 2 of the License, or
   8  * (at your option) any later version.
   9  *
  10  * This program is distributed in the hope that it will be useful, but
  11  * WITHOUT ANY WARRANTY; without even the implied warranty of
  12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13  * General Public License for more details.
  14  *
  15  * You should have received a copy of the GNU General Public License
  16  * along with this program; if not, write to the Free Software
  17  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
  18  * USA
  19  */
  20
  21 /** @file lib/bits.c
  22  * @brief Bit operations
  23  */
  24
  25 #include <config.h>
  26 #include "types.h"
  27
  28 #include <math.h>
  29
  30 #include "bits.h"
  31
  32 #if !HAVE_FLS
  33 /** @brief Compute index of leftmost 1 bit
  34  * @param n Integer
  35  * @return Index of leftmost 1 bit or -1
  36  *
  37  * For positive @p n we return the index of the leftmost bit of @p n.  For
  38  * instance @c leftmost_bit(1) returns 0, @c leftmost_bit(15) returns 3, etc.
  39  *
  40  * If @p n is zero then -1 is returned.
  41  */
  42 int leftmost_bit(uint32_t n) {
  43   /* See e.g. Hacker's Delight s5-3 (p81) for where the idea comes from.
  44    * Warren is computing the number of leading zeroes, but that's not quite
  45    * what I wanted.  Also this version should be more portable than his, which
  46    * inspects the bytes of the floating point number directly.
  47    */
  48   int x;
  49   frexp((double)n, &x);
  50   /* This gives: n = m * 2^x, where 0.5 <= m < 1 and x is an integer.
  51    *
  52    * If we take log2 of either side then we have:
  53    *    log2(n) = x + log2 m
  54    *
  55    * We know that 0.5 <= m < 1 => -1 <= log2 m < 0.  So we floor either side:
  56    *
  57    *    floor(log2(n)) = x - 1
  58    *
  59    * What is floor(log2(n))?  Well, consider that:
  60    *
  61    *    2^k <= z < 2^(k+1)  =>  floor(log2(z)) = k.
  62    *
  63    * But 2^k <= z < 2^(k+1) is the same as saying that the leftmost bit of z is
  64    * bit k.
  65    *
  66    *
  67    * Warren adds 0.5 first, to deal with the case when n=0.  However frexp()
  68    * guarantees to return x=0 when n=0, so we get the right answer without that
  69    * step.
  70    */
  71   return x - 1;
  72 }
  73 #endif
  74
  75 /*
  76 Local Variables:
  77 c-basic-offset:2
  78 comment-column:40
  79 fill-column:79
  80 indent-tabs-mode:nil
  81 End:
  82 */