/////////////////////////////////////////////////////////////////////////////
//
// Image class containing statistical function
//
// Copyright © 1998
//
// IMM, Department of Mathematical Modelling
// Technical University of Denmark, Building 321
// DK-2800 Lyngby, Denmark
// http://www.imm.dtu.dk
//
// author: Rune Fisker 
//
// Disclaimer:
//
// No guarantees of performance accompany this software,
// nor is any responsibility assumed on the part of the author(s).
// The software has been tested extensively and every effort has been
// made to insure its reliability.
//
// This software is provided by IMM and the contributor(s) ``as is'' and
// any express or implied warranties, including, but not limited to, the
// implied warranties of merchantability and fitness for a particular purpose
// are disclaimed.  In no event shall IMM or the contributor(s) be liable
// for any direct, indirect, incidental, special, exemplary, or consequential
// damages (including, but not limited to, procurement of substitute goods
// or services; loss of use, data, or profits; or business interruption)
// however caused and on any theory of liability, whether in contract, strict
// liability, or tort (including negligence or otherwise) arising in any way
// out of the use of this software, even if advised of the possibility of
// such damage.
//
// This software is partly based on the Microsoft Vision Software Developers Kit VisSDK
//
/////////////////////////////////////////////////////////////////////////////

//#include "stdafx.h"

/////////////////////////////////////////////////////////////////////////////
//	compute mean of pixels (only gray scale). 
//  Different from rgb, because of different return type (and missing rgb-type conversions)
//
//  dMean: mean value
//
//	author:	Rune Fisker, 2/9-1998
template <class TPixel>
void CDImageStatistic<TPixel>::Mean(double& dMean) 
{

	assert(NBands() == 1);
	assert(!(PixFmt() & evispixfmtRGBA));

	dMean = 0;
	for (int r = Top();r < Bottom(); r++)
	{
		TPixel *ptPixel = RowPointer(r);
		for (int c = Left();c < Right();c++)
		{
			dMean += ptPixel[c];
		}
	}

	dMean /= NPixels();
}


/////////////////////////////////////////////////////////////////////////////
//	compute mean of pixels (only rgb). 
//  Different from gray scale, because of different return type (and missing rgb-type conversions)
//
//  rgbMean: mean value
//
//	author:	Rune Fisker, 2/9-1998
template <class TPixel>
void CDImageStatistic<TPixel>::Mean(CVisRGBADoublePixel& rgbMean) 
{

	assert(NBands() == 1);

	assert(PixFmt() & evispixfmtRGBA);

	rgbMean = 0;
	rgbMean.SetA(0);
	
	for (int r = Top();r < Bottom(); r++)
	{
		TPixel *ptPixel = RowPointer(r);
		for (int c = Left();c < Right();c++)
		{			
			rgbMean.SetR(rgbMean.R() + ptPixel[c].R());
			rgbMean.SetG(rgbMean.G() + ptPixel[c].G());
			rgbMean.SetB(rgbMean.B() + ptPixel[c].B());
			rgbMean.SetA(rgbMean.A() + ptPixel[c].A());
		}
	}

	rgbMean /= NPixels();
}

/////////////////////////////////////////////////////////////////////////////
//	returns max of pixels (only gray scale). 
//
//	author:	Rune Fisker, 2/9-1998
//	Modified: Rune Fisker, 12/1-1999
template <class TPixel>
TPixel CDImageStatistic<TPixel>::Max(void) 
{

	assert(NBands() == 1);

	TPixel tMax = Pixel(Left(),Top());
	for (int r = Top();r < Bottom(); r++)
	{
		TPixel *ptPixel = RowPointer(r);
		for (int c = Left();c < Right();c++)
		{
			if (tMax < ptPixel[c])
				tMax = ptPixel[c];
		}
	}

	return tMax;
}

/////////////////////////////////////////////////////////////////////////////
//	returns min of pixels (only gray scale). 
//
//	author:	Rune Fisker, 2/9-1998
//	Modified: Rune Fisker, 12/1-1999
template <class TPixel>
TPixel CDImageStatistic<TPixel>::Min() 
{

	assert(NBands() == 1);

	TPixel tMin = Pixel(Left(),Top());
	for (int r = Top();r < Bottom(); r++)
	{
		TPixel *ptPixel = RowPointer(r);
		for (int c = Left();c < Right();c++)
		{
			if (ptPixel[c] < tMin)
				tMin = ptPixel[c];
		}
	}

	return tMin;
}


/////////////////////////////////////////////////////////////////////////////
//	compute the central (non-spatial) moment of the pixels.
//  Different from rgb, because of different return type (and missing rgb-type conversions)
//
//	input:
//		nOrder: order if the momemt
//
// output:
//		dMoment: central moment of order nOrder
//
//	author:	Rune Fisker, 19/4-1999
template <class TPixel>
void CDImageStatistic<TPixel>::CentralMomentPixel(double& dMoment, int nOrder)
{
	assert(NBands() == 1);
	assert(!(PixFmt() & evispixfmtRGBA));

	double dMean;
	Mean(dMean);

	dMoment = 0;
	for (int r = Top();r < Bottom(); r++)
	{
		TPixel *ptPixel = RowPointer(r);
		for (int c = Left();c < Right();c++)
		{
			dMoment += pow(ptPixel[c] - dMean,nOrder);
		}
	}
}

/////////////////////////////////////////////////////////////////////////////
//	compute variance of the pixels (only gray scale). 
//  Different from rgb, because of different return type (and missing rgb-type conversions)
//
//  dStd: standard deviation
//
//	author:	Rune Fisker, 19/4-1999
//	modified:	Rune Fisker, 15/7-1999 
template <class TPixel>
void CDImageStatistic<TPixel>::Var(double& dVar) 
{
	CentralMomentPixel(dVar,2);
	dVar /= NPixels();
}

/////////////////////////////////////////////////////////////////////////////
//	compute standard deviation of the pixels (only gray scale). 
//  Different from rgb, because of different return type (and missing rgb-type conversions)
//
//  dStd: standard deviation
//
//	author:	Rune Fisker, 19/4-1999
//	DIVAUPDATE
template <class TPixel>
void CDImageStatistic<TPixel>::Std(double& dStd) 
{
	Var(dStd);
	dStd = sqrt(dStd);
}

/////////////////////////////////////////////////////////////////////////////
//	compute skewnees of the pixels (only gray scale). 
//  Different from rgb, because of different return type (and missing rgb-type conversions)
//
//  dStd: standard deviation
//
//	author:	Rune Fisker, 19/4-1999
template <class TPixel>
void CDImageStatistic<TPixel>::Skewness(double& dSkewness) 
{
	double dStd;
	Std(dStd);

	CentralMoment(dSkewness);
	dSkewness /= (pow(dStd,3)*NPixels());
}