// GaussMatrix.cpp: implementation of the CGaussMatrix class.
//
//////////////////////////////////////////////////////////////////////

#include "GaussMatrix.h"
#include <stdio.h>
#include <math.h>
#define M_PI 3.1415

//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

CGaussMatrix::CGaussMatrix()
{
	mu = 0.0;
	sigma = 1.0;
	Elements = new double *[1];
	Elements[0] = new double [1];
	SizeX = 1;
	SizeY = 1;
	Elements[0][0] = 1.0;
	printf("CGaussMatrix class (new Gaussmatrix): succesfull...\n");	
}

CGaussMatrix::~CGaussMatrix()
{
	
}

CGaussMatrix::CGaussMatrix(int size)
{

	int i,j;
	double tav, norm = 0.0, const1;
	if (size % 2 != 1) size ++;
	mu = 0.0;
	sigma = 1.0;
	Elements = new double *[size];
	for (i = 0; i < size; i ++)
		Elements[i] = new double [size];
	SizeX = size;
	SizeY = size;

	
	//compute the Gauss 2 dimensional matrix

    const1 = 1 / (sqrt(sigma*2)*M_PI);

    for (i=0;i<size;i++){
        for (j=0;j<size;j++){
            tav = sqrt((-((size-1)/2)+i)*(-((size-1)/2)+i)+(-((size-1)/2)+j)*(-((size-1)/2)+j));
            Elements[i][j] = const1 * exp( - ( ( (tav - mu)*(tav - mu) ) / (2*sigma*sigma) ));
        }
    }
    i=0;

    for (i=0;i<size;i++){
        for (j=0;j<size;j++){
            norm += Elements[i][j];
        }
    }
    for (i=0;i<size;i++){
        for (j=0;j<size;j++){
            Elements[i][j] *= (1/norm);
        }
    }

	printf("CGaussMatrix class (new Gaussmatrix %d x %d): succesfull...\n", size, size);	
}

CGaussMatrix::CGaussMatrix(int size, double muf, double sigmaf)
{

	int i,j;
	double tav, norm = 0.0, const1;
	if (size % 2 != 1) size ++;
	mu = muf;
	sigma = sigmaf;
	Elements = new double *[size];
	for (i = 0; i < size; i ++)
		Elements[i] = new double [size];
	SizeX = size;
	SizeY = size;

	
	//compute the Gauss 2 dimensional matrix

    const1 = 1 / (sqrt(2*M_PI)*sigma);

    for (i=0;i<size;i++){
        for (j=0;j<size;j++){
            tav = sqrt((-((size-1)/2)+i)*(-((size-1)/2)+i)+(-((size-1)/2)+j)*(-((size-1)/2)+j));
			Elements[i][j] = const1 * exp( - ( ( (tav - mu)*(tav - mu) ) / (2*sigma*sigma) ));
        }
    }
    i=0;

    for (i=0;i<size;i++){
        for (j=0;j<size;j++){
            norm += Elements[i][j];
        }
    }
    for (i=0;i<size;i++){
        for (j=0;j<size;j++){
            Elements[i][j] *= (1/norm);
        }
    }

	printf("CGaussMatrix class (new Gaussmatrix %d x %d): succesfull...\n", size, size);	

}

double CGaussMatrix::GetXY(int x, int y)
{
	int i, j;
	i = x + SizeX/2;
	j = y + SizeY/2;
	if (i>=0 && j>= 0 && i < SizeX && j < SizeY) {
		return Elements[i][j];
	}else{
		return 0.0;
	}	
}