Thinning.cpp

/*
This file is part of the BIAS library (Basic ImageAlgorithmS).

Copyright (C) 2003-2009    (see file CONTACT for details)
  Multimediale Systeme der Informationsverarbeitung
  Institut fuer Informatik
  Christian-Albrechts-Universitaet Kiel


BIAS is free software; you can redistribute it and/or modify
it under the terms of the GNU Lesser General Public License as published by
the Free Software Foundation; either version 2.1 of the License, or
(at your option) any later version.

BIAS is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License
along with BIAS; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
*/


#include "Thinning.hh"

#ifdef BIAS_DEBUG
#include "Base/Image/ImageIO.hh"
#endif

#include <Base/Common/BIASpragma.hh>

using namespace BIAS;
using namespace std;

template <class InputStorageType, class OutputStorageType>
Thinning<InputStorageType, OutputStorageType>::
Thinning()
  : FilterNToN<InputStorageType, OutputStorageType>() 
{ 
}

template <class InputStorageType, class OutputStorageType>
Thinning<InputStorageType, OutputStorageType>::
Thinning(const Thinning<InputStorageType, OutputStorageType>& other)
    : FilterNToN<InputStorageType, OutputStorageType>(other) 
{
}

template <class InputStorageType, class OutputStorageType>
Thinning<InputStorageType, OutputStorageType>::
~Thinning()
{
}

template <class InputStorageType, class OutputStorageType>
int Thinning<InputStorageType, OutputStorageType>::
FilterInt(const Image<InputStorageType>& src, Image<OutputStorageType>& dst)
{
  BIASERR("Unimplemented interface, calling generic Filter");
  return Filter(src, dst);
}

template <class InputStorageType, class OutputStorageType>
int Thinning<InputStorageType, OutputStorageType>::
FilterFloat(const Image<InputStorageType>& src, Image<OutputStorageType>& dst)
{
  BIASERR("Unimplemented interface, calling generic Filter");
  return Filter(src, dst);
}


template <class InputStorageType,class OutputStorageType>
int Thinning<InputStorageType,OutputStorageType>::
Filter(const Image<InputStorageType> &src, Image<OutputStorageType> &dst) 
{
  BIASCDOUT(D_FILTERBASE_CALLSTACK,"Thinning::Rosenfeld_\n");
  const InputStorageType **srcstart=src.GetImageDataArray();
  OutputStorageType **dststart=dst.GetImageDataArray();
  int img_width=src.GetWidth();
  int img_height=src.GetHeight();
  int channelcount=src.GetChannelCount();

  std::vector < std::vector < bool > > Econnected_;
  std::vector < std::vector < bool > > Eendpoint_;
  std::vector < std::vector < bool > > Eisolated_;
  std::vector < std::vector < bool > > Esimple_;

  std::vector < std::vector < bool > > northBorder_;
  std::vector < std::vector < bool > > southBorder_;
  std::vector < std::vector < bool > > eastBorder_;
  std::vector < std::vector < bool > > westBorder_;

  //122 patterns for all possible not-simple variants for our CELLS
  int not_simple[] = {144,136,132,80,72,68,66,65,36,34,
                      33,20,18,17,9,5,208,200,196,164,
                      148,152,146,145,140,137,133,100,98,97,
                      84,82,81,88,76,74,73,70,69,67,
                      52,50,49,41,38,37,35,25,22,21,
                      19,13,228,216,212,210,209,204,201,197,
                      116,114,113,105,102,101,99,180,165,156,
                      153,150,149,147,141,92,90,89,86,85,
                      83,78,77,75,71,57,54,53,51,45,
                      39,29,23,220,217,214,213,211,205,241,
                      229,181,157,151,125,119,95,221,215,61,
                      55,121,118,117,115,109,103,94,93,91,
                      87,79};

  //Reserve place in these arrays
  for (int i = 0; i < img_height; i++)
    {
      std::vector< bool > just;
      just.clear();

      for (int j = 0; j < img_width; j++)
        {
          //bEdges[i].reserve(size);
          just.push_back(false);
        }

      Econnected_.push_back(just);
      Eendpoint_.push_back(just);
      Eisolated_.push_back(just);
      Esimple_.push_back(just);

      northBorder_.push_back(just);
      southBorder_.push_back(just);
      eastBorder_.push_back(just);
      westBorder_.push_back(just);
    }

  const InputStorageType **srcida = NULL;
  OutputStorageType **dstida = NULL;

  if(channelcount == 1)
    {
      srcida = srcstart;
      dstida = dststart;

      //copy image
      for(int x = 0; x < img_height; x++)
        {
          for(int y = 0; y < img_width; y++) 
            {
              dstida[x][y] = (OutputStorageType)srcida[x][y];
            }
        }
      
      //Do RosenFeld 4 times - for each border one time
      for(int bord = 0; bord < 4; bord++)
        {
          
          for(int x = 0; x < img_height; x++)
            {
              for(int y = 0; y < img_width; y++) 
                {                  
                  int sum8 = 0;
                  int sum8BIG = 0;
                  
                  if((int)dstida[x][y] != 0)
                    {
                      if(x > 0)
                        {
                          // x0x
                          // xcx
                          // xxx
                          if((int)dstida[x-1][y] > 0){
                            sum8++; sum8BIG += 64;}
                          else{
                            northBorder_[x][y] = true;}
                          
                          // 0xx
                          // xcx
                          // xxx
                          if((y > 0) && ((int)dstida[x-1][y-1] > 0)){
                            sum8++; sum8BIG += 128;}
                          
                          // xx0
                          // xcx
                          // xxx
                          if((y < img_width - 1) && ((int)dstida[x-1][y+1] > 0)){
                            sum8++; sum8BIG += 32;}
                        }
                      if(x < img_height - 1)
                        {
                          // xxx
                          // xcx
                          // x0x
                          if((int)dstida[x+1][y] > 0){
                            sum8++; sum8BIG += 4;}
                          else{
                            southBorder_[x][y] = true;}
                          
                          // xxx
                          // xcx
                          // 0xx
                          if((y > 0) && ((int)dstida[x+1][y-1] > 0)){
                            sum8++; sum8BIG += 2;}
                          
                          // xxx
                          // xcx
                          // xx0
                          if((y < img_width - 1) && ((int)dstida[x+1][y+1] > 0)){
                            sum8++; sum8BIG += 8;}
                        }
                      
                      if(y > 0)
                        {
                          // xxx
                          // 0cx
                          // xxx
                          if((int)dstida[x][y-1] > 0){
                            sum8++; sum8BIG += 1;}
                          else{
                            westBorder_[x][y] = true;}
                        }
                      
                      if(y < img_width - 1)
                        {
                          // xxx
                          // xc0
                          // xxx
                          if((int)dstida[x][y+1] > 0){
                            sum8++; sum8BIG += 16;}
                          else{
                            eastBorder_[x][y] = true;}
                        }
                    }
                  
                  if(sum8 > 0)
                    {
                      if(sum8 == 1)
                        {
                          Eendpoint_[x][y] = true;
                        }
                      
                      Econnected_[x][y] = true;
                      Esimple_[x][y] = true;
                      for(int i = 0; i < 122; i++)
                        {
                          if(sum8BIG == not_simple[i]){Esimple_[x][y] = false;}
                        }
                    }
                  else
                    {
                      Eisolated_[x][y] = true;
                    }
                }
            }
          
          for(int x = 0; x < img_height; x++)
            {
              for(int y = 0; y < img_width; y++) 
                {
                  
                  /*cout << "Point with coords: " << y << "; " << x << " has value : " << (int)srcida[x][y] 
                    << " and (con, iso, sim, end, n, s, w, e): "
                    << Econnected_[x][y] << Eisolated_[x][y] << Esimple_[x][y] << Eendpoint_[x][y]
                    << northBorder_[x][y] << southBorder_[x][y] << westBorder_[x][y] << eastBorder_[x][y] << endl;*/
                  
                  if(Esimple_[x][y] && 
                     !Eendpoint_[x][y] && 
                     !Eisolated_[x][y] &&
                     (((bord == 0) && northBorder_[x][y]) ||
                      ((bord == 1) && southBorder_[x][y]) ||
                      ((bord == 2) && eastBorder_[x][y]) ||
                      ((bord == 3) && westBorder_[x][y])
                     ))
                    {
                      dstida[x][y]=(OutputStorageType)(0);
                    }
                }
            }
        }
    }
  return 0;
}

namespace BIAS{
#define FILTER_INSTANTIATION_CLASS Thinning
#include "Filterinst.hh"
}