mrf.cpp

Supervised Image Segmentation Using Markov Random Fields

/******************************************************************
 * Modul name : mrf.cpp
 * Author     : Csaba Gradwohl ([e-mail]) with some
 *              minor contributions from Zoltan Kato ([e-mail]).
 * Copyright  : GNU General Public License www.gnu.org/copyleft/gpl.html
 * Description:
 * Intensity-based image segmentation using a Markov random field
 * segmentation model and four different optimization algorithms:
 * Metropolis - Simulated Annealing using Metropolis dynamics
 * Gibbs      - Simulated Annealing using a Gibbs sampler
 * ICM        - Iterated Conditional Modes, a deterministic suboptimal
 *              method (depends on a good initialization).
 * MMD        - Modified Metropolis Dynamics, a pseudo-stochastic
 *              suboptimal method which is less sensitive to
 *              initialization than ICM.
 *
 * The program GUI is written in wxWindows hence the code can be
 * compiled and ran under Windows as well as under Linux/Unix.
 *
 * $Id: mrf.cpp,v 1.7 2005/02/15 21:46:16 kato Exp $
 * $Revision: 1.7 $
 * $State: Exp $
 * $Log: mrf.cpp,v $
 * Revision 1.7  2005/02/15 21:46:16  kato
 * CPU timer now doesn't count GUI overhead
 *
 * Revision 1.6  2005/02/14 15:25:38  kato
 * fixed frame size
 *
 * Revision 1.5  2004/12/13 13:35:04  kato
 * Adjusted position of version & copyright info
 *
 * Revision 1.4  2004/12/13 13:31:12  kato
 * Added copyright & version information,
 * Class parameter window is larger
 *
 * Revision 1.3  2004/12/13 13:15:02  kato
 * Bug fixes
 *
 * Revision 1.2  2004/12/13 12:07:28  kato
 * added CPU timer
 *
 * Revision 1.1  2004/12/08 13:15:06  kato
 * Initial revision
 *
 * 
 *****************************************************************/
#ifndef lint
static char rcsid_mrf_cpp[]="$Id: mrf.cpp,v 1.7 2005/02/15 21:46:16 kato Exp $";
#endif

/* wxWindows includes
 */
#include <wx/wxprec.h>
#ifndef WX_PRECOMP
    #include <wx/wx.h>
#endif
#include <wx/image.h>

#include <math.h>
#include <stdlib.h>

/* Random number generators
 */
#include "randomc.h"   // define classes for random number generators

/* Timer classes
 */
#include "CKProcessTimeCounter.h"

#define WINDOW_TITLE "MRF Image Segmentation Demo $Revision: 1.7 $"
#define VERSION      "MRF Image Segmentation Demo $Revision: 1.7 $ (Last built "\
                     __DATE__" "__TIME__") "
#define COPYRIGHT    "(c) 2004 by Csaba Gradwohl & Zoltan Kato (SZTE - Hungary)"


static wxTextCtrl *gaussians;            // output textfield for Gaussian parameters
static CKProcessTimeCounter timer("core"); // CPU timer
static bool timer_valid = FALSE;

/* Program's application class 
 */
class MyApp: public wxApp
{
  virtual bool OnInit(); // this is the main entry point
};


/* ImageOperations class: it handles all image operations such as
 * loading, saving, etc... 
 */
class ImageOperations
{
public:
  ImageOperations(wxWindow *_frame);    // constructor
  wxImage *LoadBmp(wxString bmp_name);	// loads an image from file		
  bool SaveBmp(wxString bmp_name);      // saves out_image to a given file
  bool IsOutput();			// TRUE if  out_image <> NULL
  void SetNoRegions(int n);	      	// sets the number of regions,
					// allocates/frees memory for
					// means and variances
  int GetNoRegions() { return no_regions; }
  void SetBeta(double b) { beta = b; }
  void SetT(double x) { t = x; }
  void SetT0(double t) { T0 = t; }
  void SetC(double x) { c = x; }
  void SetAlpha(double x) { alpha = x; }
  int GetK() { return K; }
  double GetT() { return T; }
  double GetE() { return E; }
  double GetTimer() { return (timer_valid? timer.GetElapsedTimeMs() : 0.0); }

  void CalculateMeanAndVariance(int region);  // computes mean and
					      // variance of the given region.
  double CalculateEnergy();                   // computes global energy
					      // based on the current
					      // lableing in data
  double LocalEnergy(int i, int j, int label);// computes the local
					      // energy at site (i,j)
					      // assuming "label" has
					      // been assigned to it.

  void Metropolis(bool mmd=false);  // executes Metropolis or MMD (if mmd=true)
  void ICM();			    // executes ICM
  void Gibbs();			    // executes Gibbs sampler

private:
  wxWindow *frame;		    // the main window
  wxImage *in_image, *out_image;    // input & output images
  int width, height;		    // width and height of the
				    // displayed image
  int no_regions;	            // number of regions for Gaussian
				    // parameter computation
  double beta;                      // strength of second order clique potential
  double t;			    // Stop criteraia threshold: stop
				    // if (deltaE < t)
  double T0;		            // Initial temperature (not used by ICM)

  double c;			    // Temperature scheduler's factor:
				    // T(n+1)=c*T(n).
  double alpha;		            // alpha value for MMD
  double *mean;			    // computed mean values and
  double *variance;		    // variances for each region
  double E;			    // current global energy
  double E_old;			    // global energy in the prvious iteration
  double T;			    // current temperature
  int K;			    // current iteration #
  int **classes;		    // this is the labeled image
  int **in_image_data;		    // Intensity values of the input image 

  void InitOutImage();
  void CreateOutput();	           // creates and draws the output
				   // image based on the current labeling
  double Singleton(int i, int j, int label); // computes singleton
					     // potential at site
					     // (i,j) having a label "label"
  double Doubleton(int i, int j, int label); // computes doubleton
					     // potential at site
					     // (i,j) having a label "label"
};


/* MyScrolledWindow class: the window used for diaplaying images
 */
class MyScrolledWindow: public wxScrolledWindow
{
public:
  MyScrolledWindow(wxWindow* parent, wxWindowID id = -1, 
		   const wxPoint& pos = wxDefaultPosition, 
		   const wxSize& size = wxDefaultSize, 
		   long style = wxHSCROLL | wxVSCROLL, 
		   const wxString& name = "scrolledWindow"): 
    wxScrolledWindow(parent, id, pos, size, style, name) 
  { bmp = NULL; }
  void SetBmp(wxImage *_bmp);     // assigns the image to the window.

protected:
  virtual void OnDraw(wxDC& dc);  // displays the image in the window

private:
  wxImage *bmp;			  // the image to be displayed
  int xDst, yDst;		  // the position of the image within
				  // the window (meaningful only when
				  // the image is smaller than the window)

  wxMemoryDC memDC;		  // memDC storing the image

  void OnLeftDown(wxMouseEvent& event);	   // Left button event handler
  void OnMouseMotion(wxMouseEvent& event); // mouse motion event handler
  DECLARE_EVENT_TABLE()
};


/* MyFrame class: the main window
 */
class MyFrame: public wxFrame
{
public:
  MyFrame(const wxString& title, const wxPoint& pos, const wxSize& size);
  ~MyFrame();
  /* returns the coordinates of the given training rectangle (or the
   * current one if region==-1)
   */
  void GetRegion(int &x, int &y, int &w, int &h, int region=-1);
  
  int GetActRegion() { return act_region;}
  void SetRegs1(int x, int y) { 
    regs[act_region*4] = x; 
    regs[act_region*4+1] = y; 
  }
  void SetRegs2(int x, int y) { 
    regs[act_region*4+2] = x; 
    regs[act_region*4+3] = y; 
  }
  MyScrolledWindow *GetInputWindow() { return input_window; }
  MyScrolledWindow *GetOutputWindow() { return output_window; }
  
  bool IsSelected(int region) { // tells whether the region has been selected 
    return (regs[region*4] || regs[region*4+1] || 
	    regs[region*4+2] || regs[region*4+3]); 
  } 
  bool AllRegionsSelected() { // tells whether all the regions has been selected
    for(int i=1; i<imageop->GetNoRegions(); ++i)
      if (!IsSelected(i)) return false; 
    return true; 
  } 
  
private:
  ImageOperations *imageop; 
  MyScrolledWindow *input_window, *output_window;    // input & output
						     // images' window
  wxButton *load_button, *save_button, *doit_button; // buttons
  wxButton *select_region_button;
  wxChoice *op_choice;		// scroll-list of optimization algorithms
  wxTextCtrl *regions;          // input field for number of classes,
  wxTextCtrl *tbeta, *tt;	// beta, threshold t,
  wxTextCtrl *tT0, *tc;		// initial temperature T0, scheduler factor c,
  wxTextCtrl *talpha;		// and MMD's alpha
  int act_region;   // the current class
  int *regs;	    // stores the training rectangles for each class.
  
  /* Event handlers
   */  
  void OnOpen(wxCommandEvent& event);         // Load
  void OnSave(wxCommandEvent& event);         // Save
  void OnDoit(wxCommandEvent& event);         // DoIt
  void OnChoice(wxCommandEvent& event);	      // optimization method selection 
  void OnRegions(wxCommandEvent& event);      // number of classes
  void OnSelectRegion(wxCommandEvent& event); // select training rectangle
  void OnPaint(wxPaintEvent& event);	      // paint handler
  DECLARE_EVENT_TABLE()
    };

enum { ID_LOAD_BUTTON, ID_SAVE_BUTTON, ID_DOIT_BUTTON, ID_CHOICE,
       ID_REGIONS, ID_SELECTREGION_BUTTON, ID_BETA, ID_T, ID_T0, ID_C,
       ID_ALPHA, ID_GAUSSIANS };

/* Event table
 */
BEGIN_EVENT_TABLE(MyFrame, wxFrame)
  EVT_BUTTON(ID_LOAD_BUTTON, MyFrame::OnOpen)
  EVT_BUTTON(ID_SAVE_BUTTON, MyFrame::OnSave)
  EVT_BUTTON(ID_DOIT_BUTTON, MyFrame::OnDoit)
  EVT_CHOICE(ID_CHOICE, MyFrame::OnChoice)
  EVT_PAINT(MyFrame::OnPaint)
  EVT_TEXT(ID_REGIONS, MyFrame::OnRegions)
  EVT_BUTTON(ID_SELECTREGION_BUTTON, MyFrame::OnSelectRegion)
END_EVENT_TABLE()

BEGIN_EVENT_TABLE(MyScrolledWindow, wxScrolledWindow)
  EVT_LEFT_DOWN(MyScrolledWindow::OnLeftDown)
  EVT_MOTION(MyScrolledWindow::OnMouseMotion)
END_EVENT_TABLE()

IMPLEMENT_APP(MyApp)

bool MyApp::OnInit()
{
  MyFrame *frame = new MyFrame( WINDOW_TITLE, 
				wxPoint(0,0), wxSize(900,620) );
  frame->Show( TRUE );
  SetTopWindow( frame );
  return TRUE;
}


void MyScrolledWindow::SetBmp(wxImage *_bmp) 
{
  xDst = yDst = 0;
  // center if image is smaller than the window
  if (_bmp!=NULL)
    {
      if (_bmp->GetWidth() < 300) xDst = (300-_bmp->GetWidth())/2;
      if (_bmp->GetHeight() < 250) yDst = (250-_bmp->GetHeight())/2;
      memDC.SelectObject(*_bmp);
    }
  bmp = _bmp; 
}

	
void MyScrolledWindow::OnDraw(wxDC& dc)
{
  if (bmp != NULL)
    {
      // determine which part of the image is visible in the window.
      int x, y;
      GetViewStart(&x, &y);
      x *= 10; y *= 10;	   // must be multiplied by ScrollUnit
      // copy the visible part into the window
      int _xDst, _yDst;
      CalcUnscrolledPosition(xDst, yDst, &_xDst, &_yDst);
      dc.Blit(_xDst, _yDst, 300, 250, &memDC, x, y);
      
      // draw the training rectangle on the input image
      MyFrame *parent = (MyFrame *)GetParent();
      if (parent->GetInputWindow() == this)
	{
	  int x1, y1, w, h;
	  parent->GetRegion(x1, y1, w, h);
	  if (x1!=0 || y1!=0 || w!=0 || h!=0)
	    {
	      wxPen pen(*wxRED_PEN);
	      wxBrush brush(*wxTRANSPARENT_BRUSH);
	      dc.SetPen(pen);
	      dc.SetBrush(brush);
	      dc.DrawRectangle(x1+xDst, y1+yDst, w, h);
	    }
	}
    }
}


/* Left mouse button event handler
 */
void MyScrolledWindow::OnLeftDown(wxMouseEvent& event)
{
  // comvert window coordinates to image coordinates
  int x, y;
  GetViewStart(&x, &y);
  x *= 10; y *= 10;	 // must be multiplied by ScrollUnit

  MyFrame *frame = (MyFrame *)GetParent();
  if (frame->GetActRegion() != -1)	// in this case regs != NULL
    {
      if (event.m_x >= xDst && event.m_x < xDst+bmp->GetWidth() &&
	  event.m_y >= yDst && event.m_y < yDst+bmp->GetHeight())
	frame->SetRegs1(event.m_x+x-xDst, event.m_y+y-yDst); // scroll added
	}
}


/* Mouse motion event handler
 */
void MyScrolledWindow::OnMouseMotion(wxMouseEvent& event)
{
  if (event.LeftIsDown())
    {
      // comvert window coordinates to image coordinates
      int x, y;
      GetViewStart(&x, &y);
      x *= 10; y *= 10;		// must be multiplied by ScrollUnit
      MyFrame *frame = (MyFrame *)GetParent();
      if (frame->GetActRegion() != -1)		// in this case regs != NULL
	{
	  if (event.m_x >= xDst && event.m_x < xDst+bmp->GetWidth() &&
	      event.m_y >= yDst && event.m_y < yDst+bmp->GetHeight())
	    frame->SetRegs2(event.m_x+x-xDst, event.m_y+y-yDst); // scroll added
	  Refresh();
		}
	}
}


/*********************************************************************
/* Functions of MyFrame class
/********************************************************************/
void MyFrame::OnPaint(wxPaintEvent& event)
{
  wxPaintDC pDC(this);

  wxString str;
  str.Printf("Number of classes:");
  pDC.DrawText(str, 20, 325);
  
  str.Printf("