inline double
sqr( double x)
{
  return x * x;
}

unsigned
bit_count( unsigned char *buf, unsigned width)
{
  int bytes = width / 8;
  int bits  = width % 8;
  int count = 0;
  int i;

  for (i = 0; i < bytes; ++i) {
    count += lut_bitcount[ buf[ i]];
  }
  for (i = 0; i < bits; ++i) {
    count += (buf[ bytes] >> (7-i)) & 0x01;
  }
  return count;
}


double
Histogram::construct( int skip, const Raster &raster)
{
  int    i;
  double sum;
  int    n = raster.height() / skip;

  if (n < m_max_count) {
    delete []m_rows;
    m_rows = new unsigned[ n];
    m_max_count = n;
  }
  m_row_count = n;
  for (i = 0; i < raster.height(); i += skip) {
    m_rows[ i] = bit_count( raster[ i], raster.width());
  }

  sum = 0.0;
  for (i = 0; i < m_row_count; ++i) {
    sum += m_rows[ i];
  }
  m_mean = sum / m_row_count;

  sum = 0.0;
  for (i = 0; i < m_row_count; ++i) {
    sum += sqr( m_rows[ i] - m_mean);
  }
  m_variance = sum / m_row_count;

  m_std_dev  = sqrt( m_variance);

  return m_std_dev;
}