/*
 *
 * Copyright (C) 2012 Francisco Madrigal
 *
 * This library is free software; you can redistribute it and/or
 * modify it under the terms of the GNU Lesser General Public
 * License as published by the Free Software Foundation; either
 * version 2.1 of the License, or (at your option) any later version.
 *
 * This library 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
 * Lesser General Public License for more details.
 
 * You should have received a copy of the GNU Lesser General Public
 * License along with this library; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
 */

typedef unsigned char uchar;

__global__ void Dx(float *Imag_dst,uchar *Imag_src,int N,int Width)
{
	int idx = blockIdx.x * blockDim.x + threadIdx.x;
	
	if (idx<N){
		if (idx%Width==0){//Cero en los bordes
			Imag_dst[idx]=0.0f;
		}
		else{
			Imag_dst[idx]=(Imag_src[idx]-Imag_src[idx-1])/255.0f;
		}
	}
}

__global__ void Dy(float *Imag_dst,uchar *Imag_src,int N,int Width)
{
	int idx = blockIdx.x * blockDim.x + threadIdx.x;
	
	if (idx<N){
		if(idx<Width){//Cero en los bordes
			Imag_dst[idx]=0.0f;
		}
		else{
			Imag_dst[idx]=(Imag_src[idx]-Imag_src[idx-Width])/255.0f;
		} 
	}


}

__global__ void MG(float *ImageMG_dev,float *ImagDx_dev,float *ImagDy_dev,int N)
{
	int idx = blockIdx.x * blockDim.x + threadIdx.x;
	
	float dx,dy;
	if (idx<N){
		dx=ImagDx_dev[idx];
		dy=ImagDy_dev[idx];
		ImageMG_dev[idx]=sqrtf(dx*dx + dy*dy);
	}

}

__global__ void harris(float *ImageDst_dev,float *ImageDx_dev,float *ImageDy_dev,int N, float k){
  int idx = blockIdx.x * blockDim.x + threadIdx.x;
	float det;
  float trc;
	if (idx<N){
		// COMPLETE CODE
		/*HERE*/
  }
}

__global__ void threshold(float *ImageDst_dev, float *ImageSrc_dev, float thrhld, int N){
	int idx = blockIdx.x * blockDim.x + threadIdx.x;
	if (idx<N){
		if (ImageDst_dev[idx] < thrhld)
			ImageDst_dev[idx] = 0;
		else
			ImageDst_dev[idx] = 255;
  }
}


__global__ void harris_window(float *ImageDst_dev,float *ImageDx_dev,float *ImageDy_dev,int imageW,int imageH,int Window, float k){
	int idx = blockIdx.x * blockDim.x + threadIdx.x;
	int idy = blockIdx.y * blockDim.y + threadIdx.y;
	int index;
	float det;
	float trc;
	if (idx<imageW && idy<imageH){

		int ix_inf=idx-Window;
		int ix_sup=idx+Window;
		int iy_inf=idy-Window;
		int iy_sup=idy+Window;

		if(ix_inf<0)
			ix_inf=0;
		if(ix_sup>imageW-1)
			ix_sup=imageW-1;
		if(iy_inf<0)
			iy_inf=0;
		if(iy_sup>imageH-1)
			iy_sup=imageH-1;

		float x2;
		float y2;
		float xy;
		float sum_x2_c = 0.0f;
		float sum_y2_c = 0.0f;
		float sum_xy = 0.0f;

		for(int i=ix_inf;i<=ix_sup;i++){
			for(int j=iy_inf;j<=iy_sup;j++){
				index=imageW*j+i;
				// COMPLETE CODE
				/*HERE*/
			}
		}
		
		det = sum_x2_c*sum_y2_c - sum_xy*sum_xy;
		trc = sum_x2_c+sum_y2_c;
		ImageDst_dev[imageW*idy+idx] = det - k*trc*trc;
	}
}

extern "C"
void CUDA_Dx(float *ImageDx_dev,uchar *Image_dev,int N,int imageW,dim3 grid,dim3 threads){
	Dx<<<grid,threads>>>(ImageDx_dev,Image_dev,N,imageW);
}

extern "C"
void CUDA_Dy(float *ImageDy_dev,uchar *Image_dev,int N,int imageW,dim3 grid,dim3 threads){
	Dy<<<grid,threads>>>(ImageDy_dev,Image_dev,N,imageW);
}

extern "C"
void CUDA_MG(float *ImageMG_dev,float *ImageDx_dev,float *ImageDy_dev,int N,dim3 grid,dim3 threads){
	MG<<<grid,threads>>>(ImageMG_dev,ImageDx_dev,ImageDy_dev,N);
}

extern "C"
void CUDA_Harris(float *ImageDst_dev,float *ImageDx_dev,float *ImageDy_dev, int N, float k, dim3 grid, dim3 threads){
  harris<<<grid,threads>>>(ImageDst_dev,ImageDx_dev,ImageDy_dev,N,k);
}

extern "C"
void CUDA_Harris_window(float *ImageDst_dev,float *ImageDx_dev,float *ImageDy_dev, int imageW,int imageH,int window, float k, dim3 grid, dim3 threads){
  harris_window<<<grid,threads>>>(ImageDst_dev,ImageDx_dev,ImageDy_dev,imageW,imageH,window,k);
}
