Twinkle of Silence
 
Search this site

Recent site activity

Sample Codes of Generic Image Library

この文章について

著者:
Nishihara Hiroaki
この文章のライセンス:
この文章のライセンスはMIT Licenseとします．くわしくはhttp://sites.google.com/site/twinkleofsilence/licenseを参照ください．

作成日時:
2009年9月25日 k-means法による画素のクラスタリングに基づいた画像減色 公開
           2009年9月28日 k-means法による画素のクラスタリングに基づいた画像減色 一部修正
2009年9月28日 領域成長法を用いた画像のラベリング 公開

この文章はBoost.GILを用いたいくつかのサンプルコードを記したものです．完全なものからは程遠いところにありますが，なにかの参考になればと思います．


k-means法による画素のクラスタリングに基づいた画像減色

画像を構成する全画素にk-means法によるクラスタリングを行います．指定した数のクラスタに分配された画素群はそれぞれ属するクラスタの中心クラスタが示す階調値へと置き換えられ，結果的に画像は指定したクラスタ数と同じ色数へと減色されます．使い方は下記を参考にしてください．

// main.cpp

#include <iostream>
#include <boost/gil/gil_all.hpp>
#include <boost/gil/extension/io/png_io.hpp>
#include "clustering.hpp"

int main (int argc, char * const argv[]) {
  using namespace boost::gil;
  rgb8_image_t input;
  png_read_image("sample.png", input);
  kmeans<rgb8_view_t>(view(input), 20, 16, view(input));
  png_write_view("output.png", view(input));
  return 0;
}

clustering.hppのコードは次のようになります．現状では，channelの扱いなどgenericでない部分が多々残っています．ちなみに，my_random.hはコード中のコメントに記した範囲での乱数を返す関数を提供するクラスを記述しています．ここには載せていませんが，実際にはBoost.randomを簡単にラップしたものです．

// clustering.hpp

#ifndef __CLUSTERING_HPP__
#define __CLUSTERING_HPP__

#include <iostream>
#include <vector>
#include <boost/gil/gil_all.hpp>
#include "my_random.h"

namespace boost {
namespace gil {


template <typename PIXEL>
const double calculate_distance(const PIXEL &lhs, const PIXEL &rhs) {
    double tmp = 0.f;
    for(int i=0; i<num_channels<PIXEL>::value; ++i) {
    double lhs_value = lhs[i];
    double rhs_value = rhs[i];
        tmp += (lhs_value - rhs_value)*(lhs_value - rhs_value);
  }
    return sqrt(tmp);
}

template <typename PIXEL>
const PIXEL calculate_average(const std::vector<const PIXEL *> &cluster) {
  if (cluster.size()==0)
    return PIXEL(0.f);
    std::vector<double> sums;
    for(int i=0; i<num_channels<PIXEL>::value; ++i)
        sums.push_back(0.f);    
    for(typename std::vector<const PIXEL *>::const_iterator iter=cluster.begin(); iter!=cluster.end(); ++iter)
        for(int i=0; i<num_channels<PIXEL>::value; ++i)
            sums[i] += (**iter)[i];
    PIXEL dst;
  typedef typename channel_type<PIXEL>::type channel_t;
    for(int i=0; i<num_channels<PIXEL>::value; ++i)
        dst[i]= static_cast<channel_t>(sums[i]/cluster.size());    
    return dst;
}

template <typename PIXEL>
const int get_nearest_center(const PIXEL &pixel,
                             const std::vector<PIXEL> &centers) {
  int center_number, min_center_number;
  // first time
  min_center_number = center_number = 0;
  typename std::vector<PIXEL>::const_iterator center_iter = centers.begin();
  double min_distance = calculate_distance(*center_iter, pixel);
  ++center_iter, ++center_number;
  // second time or later
  for(; center_iter!=centers.end(); ++center_iter) {
    double tmp_distance = calculate_distance(*center_iter, pixel);
    if(tmp_distance < min_distance) {
      min_distance = tmp_distance;
      min_center_number = center_number;
    }
    ++center_number;
  }
  return min_center_number;
}

template <typename VIEW>
int init_centers_and_clusters(const typename VIEW::const_t &src,
                              std::vector<typename VIEW::value_type> &centers,
                              std::vector<std::vector<const typename VIEW::value_type *> > &clusters,
                              const int &cluster_number) {
  typedef typename VIEW::value_type PIXEL;
  typedef typename channel_type<PIXEL>::type channel_t;
 
  for (int i=0; i<cluster_number; ++i) {
    PIXEL seed;
    if (is_same<channel_t, bits32f>()) {
      RandomReal rReal(0.f, 1.f); // [0.f, 1.f]
      for (int c=0; c<num_channels<PIXEL>::value; ++c) {
        seed[c] = rReal.rand();
      }
    } else {
      channel_t min = channel_traits<channel_t>::min_value();
      channel_t max = channel_traits<channel_t>::max_value();
      RandomInt rInt(min, max); // [min, max]
      for (int c=0; c<num_channels<PIXEL>::value; ++c) {
        seed[c] = rInt.rand();
      }
    }
    centers.push_back(seed);
    clusters.push_back(std::vector<const PIXEL *>());
  }
 
  return 0;
}

template <typename VIEW>
int sort_pixels_into_clusters(const typename VIEW::const_t &src,
                              std::vector<typename VIEW::value_type> &centers,
                              std::vector<std::vector<const typename VIEW::value_type *> > &clusters) {
  typedef typename VIEW::value_type PIXEL;

  clusters.clear();
  typename std::vector<std::vector<const PIXEL *> >::iterator cluster_iter = clusters.begin();
  for (int i=0; i<centers.size(); ++i)
    clusters.push_back(std::vector<const PIXEL *>());
  typename VIEW::const_t::iterator src_iter = src.begin();
  for(; src_iter!=src.end(); ++src_iter)
    clusters[get_nearest_center(*src_iter, centers)].push_back(&(*src_iter));
  return 0;
}

template <typename PIXEL>
int update_centers(std::vector<PIXEL> &centers,
                   std::vector<std::vector<const PIXEL *> > &clusters) {
  if (centers.size()!=clusters.size())
    return -1;

  typename std::vector<PIXEL>::iterator center_iter = centers.begin();
  typename std::vector<std::vector<const PIXEL *> >::iterator cluster_iter = clusters.begin();
  for(; cluster_iter!=clusters.end(); ++cluster_iter){
    *center_iter = calculate_average(*cluster_iter);
    ++center_iter;
  }
  return 0;
}

template <typename VIEW>
int output_result(const typename VIEW::const_t &src,
                  std::vector<typename VIEW::value_type> &centers,
                  const VIEW &dst) {
  typename VIEW::iterator dst_iter = dst.begin();
  typename VIEW::const_t::iterator src_iter = src.begin();
  for (; dst_iter!=dst.end(); ++dst_iter) {
    *dst_iter = centers[get_nearest_center(*src_iter, centers)];
    ++src_iter;
  }
    return 0;
}


////////////////////////////////////////////////////////////////////////////////
//
////////////////////////////////////////////////////////////////////////////////


template <typename VIEW>
int kmeans(const typename VIEW::const_t &src,
           const int &repeat_counts, const int &cluster_number,
           const VIEW &dst) {
  typedef typename VIEW::value_type PIXEL;

    std::vector<PIXEL> centers;
    std::vector<std::vector<const PIXEL *> > clusters;
  init_centers_and_clusters<VIEW>(src, centers, clusters, cluster_number);
 
  std::cout << "Image Segmentation using k-means Clustering" << std::endl;
  std::cout << "Repeat Count: " << repeat_counts << " time(s)" << std::endl;
  std::cout << "Number of Clusters: " << clusters.size() << std::endl << std::endl;
 
    for(int i=0; i<repeat_counts; ++i){
    std::cout << i+1 << " " << std::flush;
    sort_pixels_into_clusters<VIEW>(src, centers, clusters);
    update_centers<PIXEL>(centers, clusters);
    }
  output_result<VIEW>(src, centers, dst);
 
  std::cout << "Finish!" << std::endl << std::endl;
  return 0;
}

}}  // namespace boost::gil

#endif  // __CLUSTERING_HPP__


上段が原画像(sample.png)で，下段が結果画像(output.png)です．画像をクリックすると原寸大の画像が表示されます．

領域成長法を用いた画像のラベリング

ある画素を始点として再帰的に隣接画素を特定の条件で同じ領域の一部とみなしていく手法です．ここでは，k-means法による画素のクラスタリングに基づいた画像減色を実行した結果画像に対して，関数labeling()によってラベリングを行っています．色空間中において指定した数の色数になっていても，画像平面上においていくつの領域に分かれているのかは画像に依ります．領域ごとに一意のIDを振ってラベリングした結果をmap.pngとして出力しています．

// main.cpp

#include <iostream>
#include <boost/gil/gil_all.hpp>
#include <boost/gil/extension/io/png_io.hpp>
#include "clustering.hpp"
#include "region_growing.hpp"

int main (int argc, char * const argv[]) {
  using namespace boost::gil;
  rgb8_image_t input;
  png_read_image("sample.png", input);
  kmeans<rgb8_view_t>(view(input), 20, 16, view(input));
  png_write_view("output.png", view(input));
 
  gray16_image_t map(input.dimensions());
  labeling<rgb8_view_t, gray16_view_t>(view(input), view(map));
  png_write_view("map.png", view(map));
 
  return 0;
}

関数labeling()は関数region_growing()を呼び出すことで領域のラベリングを実行しています．関数region_growing()では，同一の領域とみなすか否かの判定を行う関数を関数ポインタとして引数で受け取ります．この例では，関数func()によって，注目画素が未処理かつ基準となる階調値と等しい場合に，同一の領域とみなしています．実際に用いる際には，このあたりをカスタマイズするとよいでしょう．領域成長法は再帰を用いて説明されることが多いようですが，ここではスタックを用いた記述をしています．

// region_growing.hpp

#ifndef __REGION_GROWING_HPP__
#define __REGION_GROWING_HPP__

#include <cassert>

namespace boost {
namespace gil {

namespace detail {

template <typename VIEW>
bool is_in(const VIEW &view, const point2<int> &coord) {
  if (coord.x>=0 && coord.x<view.width() && coord.y>=0 && coord.y<view.height())
    return true;
  return false;
}

} // namespace detail

template <typename VIEW, typename MAPVIEW>
long int region_growing(const typename VIEW::const_t &src,
                        const MAPVIEW &map,
                        const point2<int> &orign_coord,
                        bool (*judge_func)(const typename VIEW::value_type basis,
                                           const typename VIEW::value_type lhs,
                                           const typename MAPVIEW::value_type rhs),
                        const typename channel_type<MAPVIEW>::type &label) {
  assert(src.dimensions()==map.dimensions());
  BOOST_STATIC_ASSERT((boost::is_same<typename color_space_type<MAPVIEW>::type ,gray_t>::value));
 
  bool is_countuped = false;
  std::vector<point2<int> > coords;
  coords.push_back(orign_coord);
  typename VIEW::value_type value = *src.xy_at(orign_coord.x, orign_coord.y);
  long int counter = 0;
  while (coords.size()!=0) {
    point2<int> coord = coords.back();
    coords.pop_back();
    typename VIEW::const_t::locator src_loc = src.xy_at(coord.x, coord.y);
    typename MAPVIEW::locator map_loc = map.xy_at(coord.x, coord.y);
    if (!detail::is_in(src, coord) ||
        !(*judge_func)(value, *src_loc, *map_loc) ) {
      continue;
    }
    if (!is_countuped)
      is_countuped = true;
    *map_loc = label;
    ++counter;
    coords.push_back(coord+point2<int>(0,-1));
    coords.push_back(coord+point2<int>(0,1));
    coords.push_back(coord+point2<int>(-1,0));
    coords.push_back(coord+point2<int>(1,0));
  }
  return counter;
}

///////////////////////////////////////////////////////////////////////////////////////////////////
//
///////////////////////////////////////////////////////////////////////////////////////////////////

template <typename VIEW, typename MAPVIEW>
bool func(const typename VIEW::value_type basis,
                const typename VIEW::value_type lhs,
                const typename MAPVIEW::value_type map) {
  typedef typename channel_type<MAPVIEW>::type map_channel_t;

  bool is_ok = false;
  if (map!=channel_traits<map_channel_t>::max_value()) // !=UNLABELED
    return is_ok; // false
  if (basis==lhs)
    is_ok = true;
  return is_ok;
}

template <typename VIEW, typename MAPVIEW>
uint32_t labeling(const typename VIEW::const_t &src,
                  const MAPVIEW &map) {
  typedef typename channel_type<MAPVIEW>::type map_channel_t;

  const map_channel_t UNLABELED = channel_traits<map_channel_t>::max_value();
  fill_pixels(map, UNLABELED);
 
  map_channel_t count = 0;
  typename VIEW::const_t::xy_locator src_loc = src.xy_at(0,0);
  typename MAPVIEW::xy_locator map_loc = map.xy_at(0,0);
  for (int y=0; y<src.height(); ++y) {
    for (int x=0; x<src.width(); ++x) {
      long int tmp = region_growing<VIEW, MAPVIEW>(src, map,
                                                   point2<int>(x,y),
                                                   func<VIEW, MAPVIEW>,
                                                   count);
      if (tmp>0) {
        // std::cout << "no." << count << " : " << tmp << " pixel(s)." << std::endl;
        ++count;
      }
      ++src_loc.x();
      ++map_loc.x();
    }
    src_loc+=point2<std::ptrdiff_t>(-src.width(),1);
    map_loc+=point2<std::ptrdiff_t>(-map.width(),1);
  }
 
  // for visualization...
  map_channel_t range = UNLABELED / count;
  typename MAPVIEW::iterator iter = map.begin();
  for (; iter!=map.end(); ++iter)
    *iter = *iter * range;
 
  return count;
}

}}  // namespace boost::gil

#endif  // __REGION_GROWING_HPP__


上段が原画像(sample.png)で，中段が減色画像(output.png)で，下段がラベリング結果を画像にしたもの(map.png)です．画像をクリックすると原寸大の画像が表示されます．