   cutCreateTimer( &timer );
   cutStartTimer( timer );  /// Start timer
 
+  /// The Gaussian Kernel Width must be odd
   if (maxKernelWidth < 1) maxKernelWidth = 1;
   if (maxKernelWidth%2 == 0) maxKernelWidth--;
   short halfkernelsize = maxKernelWidth >> 1;
 
-///calculate gaussian mask
   float *cudaGaussKernel;
   cudaMalloc((void**)&cudaGaussKernel,kernelSize);
 
+  /// Calculate gaussian kernel
   calculateGaussianKernel<<<1,maxKernelWidth,kernelSize>>>(cudaGaussKernel, gaussianVariance, maxKernelWidth);
 
-///use calculated mask to gaussian filter
-
-  ///allocate temporary memory to image data
+  /// Allocate temporary memory to image data
   float *d_tmp;
   cudaMalloc((void**) &d_tmp, size.z*sizeof(float));
 
-  /// bind a texture to the CUDA array
+  /// Bind a texture to the CUDA array
   cudaBindTexture (NULL, gaussTexRef, cudaGaussKernel);
   CUT_CHECK_ERROR("Texture bind failed");
 
-  /// host side settable texture attributes
+  /// Host side settable texture attributes
   gaussTexRef.normalized = false;
   gaussTexRef.filterMode = cudaFilterModePoint;
 
-  /// bind a texture to the CUDA array
+  /// Bind a texture to the CUDA array
   cudaBindTexture (NULL, texRef, d_img);
   CUT_CHECK_ERROR("Texture bind failed");
 
-  /// host side settable texture attributes
+  /// Host side settable texture attributes
   texRef.normalized = false;
   texRef.filterMode = cudaFilterModePoint;
 
   kernel_1DConvolutionH_texture<<<DimGrid,DimBlock,kernelSize>>>(d_tmp,size,halfkernelsize);
 
-  ///bind temporary data texture
+  /// Bind temporary data texture
   cudaUnbindTexture(texRef);
   cudaBindTexture (NULL ,texRef, d_tmp);
 
   kernel_1DConvolutionV_texture<<<DimGrid,DimBlock,kernelSize>>>(d_img,size,halfkernelsize);
 
-  ///free allocated memory
+  /// Free allocated memory
   cudaFree(d_tmp);
   cudaUnbindTexture(texRef);
   cudaFree(cudaGaussKernel);

  
  sum.x = sum.y = 0;

#pragma unroll
  for(int k=-halfkernelsize;k<(halfkernelsize+1);k++){
    sum.x += (tex1Dfetch(texRef, pixIdx + k * (((pos.x+k)>=0)*((pos.x+k)<=(size.x-1)))) * s_gauss[k+halfkernelsize]);
    sum.y += s_gauss[k+halfkernelsize];
  }


  dim3 DimGrid(blocksPerGrid);
  int nBlocks = blocksPerGrid;

  float *d_hys;
  cudaMalloc((void**) &d_hys, (size.z*sizeof(float)));
  CUT_CHECK_ERROR("Memory hysteresis image creation failed");

