#include <stdio.h>
#include <unistd.h>
#include <time.h>

__global__ void median(int* arraydevice, int length)
{
  // sharedmem points to a memoryblock shared by all threads of the block.
  // Its size in bytes is not fixed at compile time (dynamic).
  extern __shared__ int sharedmem[];

  int nthreads = blockDim.x;
  int tid = threadIdx.x;
  int i,a,b,c,tmp;

  // Copy arraydevice to sharedmem and synchronize.
  for(i=tid; i<length; i+=nthreads) sharedmem[i] = arraydevice[i];
  __syncthreads();

  // Parallel median.
  for(i=tid+1; i<length-1; i+=nthreads)
  {
    a = sharedmem[i-1];
    b = sharedmem[i];
    c = sharedmem[i+1];

    if(a > b) { tmp=a; a=b; b=tmp; }
    if(a > c) { c=a; }
    if(b > c) { b=c; }

    arraydevice[i] = b;
  }
}

int main()
{
  int *arrayhost, *arraydevice;
  int length = 10;
  int nthreads = 3;
  int size = length*sizeof(int);
  int i;

  arrayhost = (int*)malloc(size);
  cudaMalloc(&arraydevice,size);

  srand((unsigned int)clock());
  for(i=0; i<length; i++) 
    arrayhost[i] = (int)(10.0*(float)rand() / (float)RAND_MAX);

  printf("array before filter\n");
  for(i=0; i<length; i++) printf("%d ",arrayhost[i]);

  // Launch kernel with size bytes shared memory in each block.
  cudaMemcpy(arraydevice, arrayhost, size, cudaMemcpyHostToDevice);
  median<<<1, nthreads, size>>>(arraydevice,length);
  cudaDeviceSynchronize();
  cudaMemcpy(arrayhost, arraydevice, size, cudaMemcpyDeviceToHost);

  printf("\narray after filter\n");
  for(int i=0; i<length; i++) printf("%d ",arrayhost[i]);
  printf("\n");

  return 0;
}