/*****************************************************************
 * flthin.c: FBM Release 1.2 18-Apr-93 Michael Mauldin
 *
 * Copyright (C) 1993 by Michael Mauldin.  Permission is granted
 * to use this file in whole or in part for any purpose, educational,
 * recreational or commercial, provided that this copyright notice
 * is retained unchanged.  This software is available to all free of
 * charge by anonymous FTP and in the UUNET archives.
 *
 * flthin.c: 
 *	Thin an image using the algorithm from "Digital Image
 *	Processing", page 492-3.
 *
 * CONTENTS
 *	thin_fbm (input, output)
 *	FBM *input, *output;
 *
 * EDITLOG
 *	LastEditDate = Sun Apr 18 23:37:20 1993 - Michael Mauldin
 *	LastFileName = /usr0/mlm/src/fbm/flthin.c
 *
 * HISTORY
 * 18-Apr-93  Michael Mauldin (mlm) at Carnegie-Mellon University
 *	Created from flshrp.c
 *****************************************************************/

# include <stdio.h>
# include <math.h>
# include <ctype.h>
# include "fbm.h"

/****************************************************************
 * thin_fbm: thin regions
 ****************************************************************/

#ifndef lint
static char *fbmid =
"$FBM flthin.c <1.2> 18-Apr-93  (C) 1993 by Michael Mauldin, source \
code available free from MLM@CS.CMU.EDU and from UUNET archives$";
#endif

thin_fbm (input, output)
FBM *input, *output;
{ register int i, j, k;
  register unsigned char *ibm, *obm, *tbm, *tail;
  int rowlen, w, h, round=0, deleted=1;
  FBM temp;

  if (input->hdr.physbits != 8)
  { fprintf (stderr, "thin_fbm: only 8 bit inputs may be thinned\n");
    return (0);
  }

  if (input->hdr.planes != 1 || input->hdr.clrlen > 0)
  { fprintf (stderr, "thin_fbm: only greyscale images may be thinned, use clr2gray first.\n");
    return (0);
  }

  /* Allocate output */
  output->hdr = input->hdr;
  alloc_fbm (output);

  /* Allocate temp */
  temp.hdr = input->hdr;
  alloc_fbm (&temp);

  /* Cache important size values */
  w = input->hdr.cols;
  h = input->hdr.rows;
  rowlen = input->hdr.rowlen;
  
  /* 
   * Since we dont want to change the input, and we need two arrays
   * two hold the images, we copy the input to the output, first.
   * We make sure each bit is set to 0 or 1.
   */

  for (ibm = input->bm,
		obm = output->bm,
		tail = input->bm + input->hdr.plnlen;
       ibm < tail; )
  { *obm++ = *ibm++ ? 1 : 0; }
  
  /* Do each pass repeatedly until no more points are deleted */
  while (deleted)
  { register int p2, p3, p4, p5, p6, p7, p8, p9, n0, s0;

    round++;
    deleted = 0;
    
    /*
     * Step one, remove points south and east of regions,
     * output->cm is the input, and temp.cm is the output.
     */

    for (j = 1; j < h-1; j++)
    { obm = &(output->bm[j*rowlen]);
      tbm = &(temp.bm[j*rowlen]);
      
      for (i=1; i < w-1; i++)
      { if (obm[i] == 0)
        { tbm[i] = 0; }
	
	else
	{
	  /*
	   *   	 p9 | p2 | p3
	   *	----+----+----
	   *	 p8 | p1 | p4
	   *	----+----+----
	   *	 p7 | p6 | p5
	   */
      
	  p2 = obm[i-rowlen];
	  p4 = obm[i+1];
	  p6 = obm[i+rowlen];
	  p8 = obm[i-1];
	  
	  /* book conditions (c) and (d) */
	  if (p4 && p6 && (p2 || p8))
	  { tbm[i] = 1; continue; }

	  p3 = obm[i-rowlen+1];
	  p5 = obm[i+rowlen+1];
	  p7 = obm[i+rowlen-1];
	  p9 = obm[i-rowlen-1];

	  n0 = p2 + p3 + p4 + p5 + p6 + p7 + p8 + p9;
	  
	  /* Check number of neighbors */
	  if (n0 < 2 || n0 > 6)
	  { tbm[i] = 1; continue; }
	  
	  /* Count number of transitions */
	  s0 = (((p2 == 0) && (p3 == 1)) + 
	        ((p3 == 0) && (p4 == 1)) + 
	        ((p4 == 0) && (p5 == 1)) + 
	        ((p5 == 0) && (p6 == 1)) + 
	        ((p6 == 0) && (p7 == 1)) + 
	        ((p7 == 0) && (p8 == 1)) + 
	        ((p8 == 0) && (p9 == 1)) + 
	        ((p9 == 0) && (p2 == 1)));
	  
	  if (s0 != 1)
	  { tbm[i] = 1; continue; }
	  
	  tbm[i] = 0;
	  deleted++;
	}
      }
    }
    
    /*
     * Step two, remove points north and west of regions,
     * temp.cm is the input, and output->cm is the output.
     */
    
    for (j = 1; j < h-1; j++)
    { obm = &(output->bm[j*rowlen]);
      tbm = &(temp.bm[j*rowlen]);
      
      for (i=1; i < w-1; i++)
      { if (tbm[i] == 0)
        { obm[i] = 0; }
	
	else
	{
	  /*
	   *   	 p9 | p2 | p3
	   *	----+----+----
	   *	 p8 | p1 | p4
	   *	----+----+----
	   *	 p7 | p6 | p5
	   */
      
	  p2 = tbm[i-rowlen];
	  p4 = tbm[i+1];
	  p6 = tbm[i+rowlen];
	  p8 = tbm[i-1];
	  
	  /* book conditions (c') and (d') */
	  if (p2 && p8 && (p4 || p8))
	  { obm[i] = 1; continue; }

	  p3 = tbm[i-rowlen+1];
	  p5 = tbm[i+rowlen+1];
	  p7 = tbm[i+rowlen-1];
	  p9 = tbm[i-rowlen-1];

	  n0 = p2 + p3 + p4 + p5 + p6 + p7 + p8 + p9;
	  
	  /* Check number of neighbors */
	  if (n0 < 2 || n0 > 6)
	  { obm[i] = 1; continue; }
	  
	  /* Count number of transitions */
	  s0 = (((p2 == 0) && (p3 == 1)) + 
	        ((p3 == 0) && (p4 == 1)) + 
	        ((p4 == 0) && (p5 == 1)) + 
	        ((p5 == 0) && (p6 == 1)) + 
	        ((p6 == 0) && (p7 == 1)) + 
	        ((p7 == 0) && (p8 == 1)) + 
	        ((p8 == 0) && (p9 == 1)) + 
	        ((p9 == 0) && (p2 == 1)));
	  
	  if (s0 != 1)
	  { obm[i] = 1; continue; }
	  
	  obm[i] = 0;
	  deleted++;
	}
      }
    }
    
    fprintf (stderr, "round %3d: %5d deletions\n", round, deleted);
  }
  
  for (obm = output->bm,
		tail = output->bm + output->hdr.plnlen;
       obm < tail; 
       obm++)
  { *obm = *obm ? WHITE : BLACK; }
  
  free_fbm (&temp);
  
  return (1);
}