/*
*
* Copyright (C) 2005 John Ashley Burgoyne and Ichiro Fujinaga
*               2007 Uma Kompella and Christoph Dalitz
*
* This program is free software; you can redistribute it and/or modify it
* under the terms of the GNU General Public License as published by the Free
* Software Foundation; either version 2 of the License, or (at your option)
* any later version.
*
* This program 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 General Public License for
* more details.
* 
* You should have received a copy of the GNU General Public License along
* with this program; if not, write to the Free Software Foundation, Inc., 59
* Temple Place - Suite 330, Boston, MA 02111-1307, USA.
*/

/*
 Thresholds an image according to Gatos et al.'s method. See:

 Gatos, Basilios, Ioannis Pratikakis, and Stavros
 J. Perantonis. 2004. An adaptive binarization technique for low
 quality historical documents. *Lecture Notes in Computer
 Science* 3163: 102-113.

 *background*
 Estimated background of the image.

 *binarization*
 A preliminary binarization of the image.

 Use the default settings for the other parameters unless you know
 what you are doing.

 Uses Gatos Background as an estimated background of the image.

 Uses Niblack Thresholding as a preliminary binarization

 Own parameters: double "q", default=0.6
                 double "p1", default=0.5
                 double "p2", default=0.8
*/

// This algorithm was taken from the gamera.sf.net sourcecodes
// and adopted for the FreeImage library
//
//	Copyright (C) 2007-2008:
//	monday2000	monday2000@yandex.ru

#include "FreeImage.h"
#include "Utilities.h"

extern FIBITMAP* nib_thres(FIBITMAP* src_dib, int region_size,
						   double sensitivity, int lower_bound, int upper_bound);

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

inline BYTE Get1BitPixel(BYTE *bits, unsigned x)
{	// this function is simplified from FreeImage_GetPixelIndex
	
	return (bits[x >> 3] & (0x80 >> (x & 0x07))) != 0 ? 255 : 0;
}

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

inline void SetPixel(BYTE *bits, unsigned x, BYTE* value)
{   // this function is simplified from FreeImage_SetPixelIndex
	
	*value ? bits[x >> 3] |= (0x80 >> (x & 0x7)) : bits[x >> 3] &= (0xFF7F >> (x & 0x7));
}

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

FIBITMAP* gatos_background(FIBITMAP* src_dib, FIBITMAP* bin_dib, int region_size)
{
	unsigned width = FreeImage_GetWidth(src_dib);
	
	unsigned height = FreeImage_GetHeight(src_dib);
	
	unsigned pitch = FreeImage_GetPitch(src_dib);
	
	unsigned bpp = FreeImage_GetBPP(src_dib);
	
	FIBITMAP* dst_dib = FreeImage_Allocate(width, height, bpp);
	
	if(bpp == 8)
	{
		if(FreeImage_GetColorType(src_dib) == FIC_MINISWHITE) {
			// build an inverted greyscale palette
			RGBQUAD *dst_pal = FreeImage_GetPalette(dst_dib);
			for(int i = 0; i < 256; i++) {
				dst_pal[i].rgbRed = dst_pal[i].rgbGreen =
					dst_pal[i].rgbBlue = (BYTE)(255 - i);
			}
		} else {
			// build a greyscale palette
			RGBQUAD *dst_pal = FreeImage_GetPalette(dst_dib);
			for(int i = 0; i < 256; i++) {
				dst_pal[i].rgbRed = dst_pal[i].rgbGreen =
					dst_pal[i].rgbBlue = (BYTE)i;
			}
		}
	}
	
	unsigned bin_pitch = FreeImage_GetPitch(bin_dib);
	
	BYTE* src_bits = (BYTE*)FreeImage_GetBits(src_dib); // The image raster
	
	BYTE* dst_bits = (BYTE*)FreeImage_GetBits(dst_dib); // The image raster
	
	BYTE* bin_bits = (BYTE*)FreeImage_GetBits(bin_dib); // The image raster	
	
	BYTE* lines, *linek, *lined, *lineb, *bin_linek;
	
	unsigned x, y, i, j;
	
	int k_index;
	
	BYTE* end_row = src_bits + (height-1) * pitch;
	
	BYTE* bin_end_row = bin_bits + (height-1) * bin_pitch;
	
	int end_col = width - 1;
	
	int sum;
	
	int bin_sum;
	
	
	
	for (y = 0; y < height; y++)
	{
		lines = src_bits + y * pitch;
		
		lined = dst_bits + y * pitch;
		
		lineb = bin_bits + y * bin_pitch;
		
		for (x = 0; x < width; x++)
		{
			if(Get1BitPixel(lineb, x)) // white
				
				lined[x] = lines[x];
			
			else
			{
				
				sum = 0;
				
				bin_sum = 0;
				
				// kernel processing
				for (i = 0; i < region_size; i++)
					for (j = 0; j < region_size; j++)
					{	
						// accumulate binary pixels.
						
						bin_linek = lineb + (i-1) * bin_pitch;
						
						if (bin_linek < bin_bits) bin_linek = bin_bits;
						if (bin_linek > bin_end_row) bin_linek = bin_end_row;
						
						
						// accumulate pixels.
						
						linek = lines + (i-1) * pitch;
						
						if (linek < src_bits) linek = src_bits;
						if (linek > end_row) linek = end_row;
						
						k_index = x+j-1;
						
						if (k_index < 0) k_index = 0;
						if (k_index > end_col) k_index = end_col;
						
						
						if(Get1BitPixel(bin_linek, k_index)) // white
						{					
							bin_sum++;
							
							sum += linek[k_index];
						}					
					}
					
					lined[x] = (bin_sum > 0 ? (BYTE)((double)sum / bin_sum) : 255);
					
			}			
		}
	}
	
	// Copying the DPI...
	
	FreeImage_SetDotsPerMeterX(dst_dib, FreeImage_GetDotsPerMeterX(src_dib));
	
	FreeImage_SetDotsPerMeterY(dst_dib, FreeImage_GetDotsPerMeterY(src_dib));
	
	return dst_dib;
}

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

FIBITMAP* ProcessFilter(FIBITMAP* src_dib, FIBITMAP* bg_dib, FIBITMAP* bin_dib, double q, double p1, double p2)
{
	unsigned width = FreeImage_GetWidth(src_dib);
	
	unsigned height = FreeImage_GetHeight(src_dib);
	
	unsigned pitch = FreeImage_GetPitch(src_dib);
	
	unsigned bpp = FreeImage_GetBPP(src_dib);
	
	FIBITMAP* dst_dib = FreeImage_Allocate(width, height, 1);
	
	// Build a monochrome palette
	RGBQUAD *pal = FreeImage_GetPalette(dst_dib);
	pal[0].rgbRed = pal[0].rgbGreen = pal[0].rgbBlue = 0;
	pal[1].rgbRed = pal[1].rgbGreen = pal[1].rgbBlue = 255;
	
	unsigned bin_pitch = FreeImage_GetPitch(bin_dib);
	
	BYTE* src_bits = (BYTE*)FreeImage_GetBits(src_dib); // The image raster
	
	BYTE* dst_bits = (BYTE*)FreeImage_GetBits(dst_dib); // The image raster
	
	BYTE* bg_bits = (BYTE*)FreeImage_GetBits(bg_dib); // The image raster
	
	BYTE* bin_bits = (BYTE*)FreeImage_GetBits(bin_dib); // The image raster	
	
	BYTE* lines, *lined, *lineb, *linebg;
	
	unsigned x, y;
	
	BYTE* end_row = src_bits + (height-1) * pitch;
	
	BYTE* bin_end_row = bin_bits + (height-1) * bin_pitch;
	
	int end_col = width - 1;
	
	int sum = 0;
	
	int bin_sum = 0;
	
	
	unsigned delta_numerator = 0;
	
	for (y = 0; y < height; y++)
	{
		lines = src_bits + y * pitch;
		
		linebg = bg_bits + y * pitch;
		
		for (x = 0; x < width; x++)		
			
			delta_numerator += lines[x] - linebg[x];
	}	
	
	
	unsigned delta_denominator = 0;
	
	for (y = 0; y < height; y++)
	{
		lineb = bin_bits + y * bin_pitch;
		
		for (x = 0; x < width; x++)
			
			if(!Get1BitPixel(lineb, x)) // black
				
				delta_denominator++;
	}
	
	double delta = (double)delta_numerator / delta_denominator;
	
	
	for (y = 0; y < height; y++)
	{
		lineb = bin_bits + y * bin_pitch;
		
		linebg = bg_bits + y * pitch;
		
		for (x = 0; x < width; x++)
			
			if(Get1BitPixel(lineb, x)) // white
			{					
				bin_sum++;
				
				sum += linebg[x];
			}
	}
	
	double b = (double)sum / bin_sum;
	
	
	BYTE val;
	
	for (y = 0; y < height; y++)
	{
		lines = src_bits + y * pitch;
		
		linebg = bg_bits + y * pitch;
		
		lined = dst_bits + y * bin_pitch;
		
		for (x = 0; x < width; x++)
		{
            val = (BYTE) (linebg[x] - lines[x])				
				
				> q * delta * (((1 - p2) / (1 + exp(((-4 * (int)linebg[x]) / (b * (1 - p1))) 
				+ ((2 * (1 + p1)) / (1 - p1))))) + p2)
                ? 0 : 255;
			
			SetPixel(lined, x, &val);
		}
	}
	
	// Copying the DPI...
	
	FreeImage_SetDotsPerMeterX(dst_dib, FreeImage_GetDotsPerMeterX(src_dib));
	
	FreeImage_SetDotsPerMeterY(dst_dib, FreeImage_GetDotsPerMeterY(src_dib));
	
	return dst_dib;
}

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

/**
FreeImage error handler
@param fif Format / Plugin responsible for the error 
@param message Error message
*/
void FreeImageErrorHandler(FREE_IMAGE_FORMAT fif, const char *message) {
	printf("\n*** "); 
	printf("%s Format\n", FreeImage_GetFormatFromFIF(fif));
	printf(message);
	printf(" ***\n");
}

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

/** Generic image loader

  @param lpszPathName Pointer to the full file name
  @param flag Optional load flag constant
  @return Returns the loaded dib if successful, returns NULL otherwise
*/

FIBITMAP* GenericLoader(const char* lpszPathName, int flag)
{	
	FREE_IMAGE_FORMAT fif = FIF_UNKNOWN;
	// check the file signature and deduce its format
	// (the second argument is currently not used by FreeImage)
	
	fif = FreeImage_GetFileType(lpszPathName, 0);
	
	FIBITMAP* dib;
	
	if(fif == FIF_UNKNOWN)
	{
		// no signature ?
		// try to guess the file format from the file extension
		fif = FreeImage_GetFIFFromFilename(lpszPathName);
	}
	
	// check that the plugin has reading capabilities ...
	if((fif != FIF_UNKNOWN) && FreeImage_FIFSupportsReading(fif))
	{
		// ok, let's load the file
		dib = FreeImage_Load(fif, lpszPathName, flag);
		
		// unless a bad file format, we are done !
		if (!dib)
		{
			printf("%s%s%s\n","File \"", lpszPathName, "\" not found.");
			return NULL;
		}
	}	
	
	return dib;
}

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

int main(int argc, char *argv[]) {
	
	// call this ONLY when linking with FreeImage as a static library
#ifdef FREEIMAGE_LIB
	FreeImage_Initialise();
#endif // FREEIMAGE_LIB
	
	// initialize your own FreeImage error handler
	
	FreeImage_SetOutputMessage(FreeImageErrorHandler);
	
	if(argc != 9) {
		printf("Usage : gatos_thres <input_file> <region_size (int)> <sensitivity (double)> <lower_bound (int)> <upper_bound (int)> <q (double)> <p1 (double)> <p2 (double)>\n");
		return 0;
	}
	
	FIBITMAP *dib = GenericLoader(argv[1], 0);	
	
	int region_size = atoi(argv[2]);
	
	double sensitivity = atof(argv[3]);
	
	int lower_bound = atoi(argv[4]);
	
	int upper_bound = atoi(argv[5]);
	
	double q = atoi(argv[6]);
	
	double p1 = atoi(argv[7]);
	
	double p2 = atoi(argv[8]);
	
	if (dib)
	{		
		// bitmap is successfully loaded!
		
		if (FreeImage_GetImageType(dib) == FIT_BITMAP)
		{
			if (FreeImage_GetBPP(dib) == 8)
			{
				FIBITMAP* bin_dib = nib_thres(dib, region_size, sensitivity, lower_bound, upper_bound);				
				
				FIBITMAP* bg_dib = gatos_background(dib, bin_dib, region_size);
				
				FIBITMAP* dst_dib = ProcessFilter(dib, bg_dib, bin_dib, q, p1, p2);
				
				if (dst_dib)
				{					
					// save the filtered bitmap
					const char *output_filename = "filtered.tif";
					
					const char *bin_output_filename = "niblack_treshold.tif";
					
					const char *bg_output_filename = "gatos_background.tif";
					
					// first, check the output format from the file name or file extension
					FREE_IMAGE_FORMAT out_fif = FreeImage_GetFIFFromFilename(output_filename);
					
					if(out_fif != FIF_UNKNOWN)
					{
						// then save the file
						FreeImage_Save(out_fif, bin_dib, bin_output_filename, 0);
						
						FreeImage_Save(out_fif, bg_dib, bg_output_filename, 0);
						
						FreeImage_Save(out_fif, dst_dib, output_filename, 0);
					}
					
					// free the loaded FIBITMAP
					FreeImage_Unload(dst_dib);					
				}			
				
				FreeImage_Unload(bin_dib);
				
				FreeImage_Unload(bg_dib);
			}
			
			else
				
				printf("%s\n", "Unsupported color mode.");
		}
		
		else // non-FIT_BITMAP images are not supported.
			
			printf("%s\n", "Unsupported color mode.");
		
		FreeImage_Unload(dib);
	}	 
	
	// call this ONLY when linking with FreeImage as a static library
#ifdef FREEIMAGE_LIB
	FreeImage_DeInitialise();
#endif // FREEIMAGE_LIB
	
	return 0;
}