# include # include # include # include int main ( int argc, char *argv[] ); void filter ( int m, int n, int *r_cpu ); __global__ void median_news_filter ( int m, int n, int r1[], int r2[], int p[] ); __device__ void i4vec_sort_bubble_a ( int n, int a[] ); void ppma_read_data ( FILE *input, int xsize, int ysize, int *r, int *g, int *b ); void ppma_read_header ( FILE *input, int *xsize, int *ysize, int *rgb_max ); int ppma_write ( char *file_out_name, int xsize, int ysize, int *r, int *g, int *b ); int ppma_write_data ( FILE *file_out, int xsize, int ysize, int *r, int *g, int *b ); int ppma_write_header ( FILE *file_out, char *file_out_name, int xsize, int ysize, int rgb_max ); /******************************************************************************/ int main ( int argc, char *argv[] ) /******************************************************************************/ /* Purpose: MAIN is the main program for DENOISE. Licensing: This code is distributed under the GNU LGPL license. Modified: 27 March 2017 Author: John Burkardt */ { int *b_cpu; int *g_cpu; char input_filename[] = "balloons_noisy.ascii.ppm"; FILE *input_unit; int m; int memsize; int n; char output_filename[] = "balloons_median.ascii.ppm"; int *r_cpu; int rgb_max; printf ( "\n" ); printf ( "DENOISE\n" ); printf ( " C/CUDA version\n" ); printf ( " Remove noise from an image.\n" ); /* Open the input file. */ input_unit = fopen ( input_filename, "r" ); if ( ! input_unit ) { fprintf ( stderr, "\n" ); fprintf ( stderr, "DENOISE - Fatal error!\n" ); fprintf ( stderr, " Could not open the file \"%s\"\n", input_filename ); exit ( 1 ); } /* Read the file header to get the data size. */ ppma_read_header ( input_unit, &m, &n, &rgb_max ); printf ( "\n" ); printf ( " Number of columns = %d\n", m ); printf ( " Number of rows = %d\n", n ); printf ( " Maximum pixel intensity = %d\n", rgb_max ); /* Allocate memory for the data. */ memsize = m * n * sizeof ( int ); r_cpu = ( int * ) malloc ( memsize ); g_cpu = ( int * ) malloc ( memsize ); b_cpu = ( int * ) malloc ( memsize ); /* Read the data. */ ppma_read_data ( input_unit, m, n, r_cpu, g_cpu, b_cpu ); fclose ( input_unit ); /* Filter R, G and B. */ filter ( m, n, r_cpu ); filter ( m, n, g_cpu ); filter ( m, n, b_cpu ); /* Write the filtered image. */ ppma_write ( output_filename, m, n, r_cpu, g_cpu, b_cpu ); printf ( "\n" ); printf ( " Wrote denoised image to \"%s\".\n", output_filename ); /* Free memory. */ free ( r_cpu ); free ( g_cpu ); free ( b_cpu ); /* Terminate. */ printf ( "\n" ); printf ( "DENOISE\n" ); printf ( " Normal end of execution.\n" ); return 0; } /******************************************************************************/ void filter ( int m, int n, int *r_cpu ) /******************************************************************************/ /* Purpose: FILTER sets up the filtering process for a single color channel. Discussion: For an RGB image, call this function once for each of R, G and B. Licensing: This code is distributed under the GNU LGPL license. Modified: 04 April 2017 Author: John Burkardt Parameters: Input, int M, N, the number of rows and columns of pixels. Input/output, int R_CPU[M*N], the color data. */ { int memsize; int *p_gpu; int *r1_gpu; int *r2_gpu; memsize = m * n * sizeof ( int ); /* Allocate R1, R2 and workspace P on the GPU. */ cudaMalloc ( ( void** ) &r1_gpu, memsize ); cudaMalloc ( ( void** ) &r2_gpu, memsize ); cudaMalloc ( ( void** ) &p_gpu, 5 * memsize ); /* Copy R values to the GPU. */ cudaMemcpy ( r1_gpu, r_cpu, memsize, cudaMemcpyHostToDevice ); /* Determine the blocking and threading. */ dim3 blocks ( 1+(m-1)/16, 1+(n-1)/16, 1 ); dim3 threads ( 16, 16, 1 ); /* Filter the R, G and B values on the GPU. */ median_news_filter <<< blocks, threads >>> ( m, n, r1_gpu, r2_gpu, p_gpu ); /* Copy the filtered GPU data back to CPU. */ cudaMemcpy ( r_cpu, r2_gpu, memsize, cudaMemcpyDeviceToHost ); /* Free memory. */ cudaFree ( p_gpu ); cudaFree ( r1_gpu ); cudaFree ( r2_gpu ); return; } /******************************************************************************/ __global__ void median_news_filter ( int m, int n, int r1[], int r2[], int p[] ) /******************************************************************************/ /* Purpose: MEDIAN_NEWS_FILTER uses a median NEWS filter on a single color array. Licensing: This code is distributed under the GNU LGPL license. Modified: 27 March 2017 Author: John Burkardt Parameters: Input, int M, N, the number of rows and columns of pixels. Input, int R1[M*N], the original color data. Output, int R2[M*N], the filtered color data. Workspace, int P[5*M*N]. */ { int k; int x; int y; x = threadIdx.x + blockIdx.x * blockDim.x; y = threadIdx.y + blockIdx.y * blockDim.y; k = x + y * m; if ( k < m * n ) { if ( x == 0 || x == m - 1 || y == 0 || y == n - 1 ) { r2[k] = r1[k]; } else { p[k*5+0] = r1[k-m]; p[k*5+1] = r1[k-1]; p[k*5+2] = r1[k]; p[k*5+3] = r1[k+1]; p[k*5+4] = r1[k+m]; i4vec_sort_bubble_a ( 5, p+k*5 ); r2[k] = p[k*5+2]; } } return; } /******************************************************************************/ __device__ void i4vec_sort_bubble_a ( int n, int a[] ) /******************************************************************************/ /* Purpose: I4VEC_SORT_BUBBLE_A ascending sorts an I4VEC using bubble sort. Discussion: An I4VEC is a vector of I4's. Licensing: This code is distributed under the GNU LGPL license. Modified: 08 June 2010 Author: John Burkardt Parameters: Input, int N, length of input array. Input/output, int A[N]. On input, an unsorted array of ints. On output, A has been sorted. */ { int i; int j; int temp; for ( i = 0; i < n - 1; i++ ) { for ( j = i + 1; j < n; j++ ) { if ( a[j] < a[i] ) { temp = a[i]; a[i] = a[j]; a[j] = temp; } } } return; } /******************************************************************************/ void ppma_read_data ( FILE *input, int xsize, int ysize, int *r, int *g, int *b ) /******************************************************************************/ /* Purpose: PPMA_READ_DATA reads the data in an ASCII portable pixel map file. Licensing: This code is distributed under the GNU LGPL license. Modified: 27 May 2011 Author: John Burkardt Parameters: Input, FILE *INPUT, a pointer to the file containing the ASCII portable pixel map data. Input, int XSIZE, YSIZE, the number of rows and columns of data. Output, int *R, *G, *B, the arrays of XSIZE by YSIZE data values. */ { int i; int j; int n; int p; p = 0; for ( j = 0; j < ysize; j++ ) { for ( i = 0; i < xsize; i++ ) { n = fscanf ( input, "%d %d %d", r, g, b ); if ( n != 3 ) { fprintf ( stderr, "\n" ); fprintf ( stderr, "PPMA_READ_DATA - Fatal error!\n" ); fprintf ( stderr, " Unable to read data.\n" ); fprintf ( stderr, " FSCANF returns N = %d\n", n ); fprintf ( stderr, " Number of pixels already read is %d\n", p ); exit ( 1 ); } p = p + 1; r = r + 1; g = g + 1; b = b + 1; } } return; } /******************************************************************************/ void ppma_read_header ( FILE *input, int *xsize, int *ysize, int *rgb_max ) /******************************************************************************/ /* Purpose: PPMA_READ_HEADER reads the header of an ASCII portable pixel map file. Licensing: This code is distributed under the GNU LGPL license. Modified: 27 May 2011 Author: John Burkardt Parameters: Input, FILE *INPUT, a pointer to the file containing the ASCII portable pixel map data. Output, int *XSIZE, *YSIZE, the number of rows and columns of data. Output, int *RGB_MAX, the maximum RGB value. */ { # define LINEMAX 255 int count; char *error; char line[LINEMAX]; char *next; int step; int width; char word[LINEMAX]; step = 0; while ( 1 ) { error = fgets ( line, LINEMAX, input ); if ( !error ) { fprintf ( stderr, "\n" ); fprintf ( stderr, "PGMA_READ_HEADER - Fatal error!\n" ); fprintf ( stderr, " End of file.\n" ); exit ( 1 ); } next = line; if ( line[0] == '#' ) { continue; } if ( step == 0 ) { count = sscanf ( next, "%s%n", word, &width ); if ( count == EOF ) { continue; } next = next + width; if ( strcmp ( word, "P3" ) != 0 && strcmp ( word, "p3" ) != 0 ) { fprintf ( stderr, "\n" ); fprintf ( stderr, "PPMA_READ_HEADER - Fatal error.\n" ); fprintf ( stderr, " Bad magic number = \"%s\".\n", word ); exit ( 1 ); } step = 1; } if ( step == 1 ) { count = sscanf ( next, "%d%n", xsize, &width ); next = next + width; if ( count == EOF ) { continue; } step = 2; } if ( step == 2 ) { count = sscanf ( next, "%d%n", ysize, &width ); next = next + width; if ( count == EOF ) { continue; } step = 3; } if ( step == 3 ) { count = sscanf ( next, "%d%n", rgb_max, &width ); next = next + width; if ( count == EOF ) { continue; } break; } } return; # undef LINEMAX } /******************************************************************************/ int ppma_write ( char *file_out_name, int xsize, int ysize, int *r, int *g, int *b ) /******************************************************************************/ /* Purpose: PPMA_WRITE writes the header and data for an ASCII portable pixel map file. Example: P3 # feep.ppm 4 4 15 0 0 0 0 0 0 0 0 0 15 0 15 0 0 0 0 15 7 0 0 0 0 0 0 0 0 0 0 0 0 0 15 7 0 0 0 15 0 15 0 0 0 0 0 0 0 0 0 Licensing: This code is distributed under the GNU LGPL license. Modified: 28 February 2003 Author: John Burkardt Parameters: Input, char *FILE_OUT_NAME, the name of the file to contain the ASCII portable pixel map data. Input, int XSIZE, YSIZE, the number of rows and columns of data. Input, int *R, *G, *B, the arrays of XSIZE by YSIZE data values. Output, int PPMA_WRITE, is true, if an error was detected, or false, if the file was written. */ { int *b_index; int error; FILE *file_out; int *g_index; int i; int j; int *r_index; int rgb_max; /* Open the output file. */ file_out = fopen ( file_out_name, "wt" ); if ( !file_out ) { printf ( "\n" ); printf ( "PPMA_WRITE - Fatal error!\n" ); printf ( " Cannot open the output file \"%s\".\n", file_out_name ); return 1; } /* Compute the maximum. */ rgb_max = 0; r_index = r; g_index = g; b_index = b; for ( j = 0; j < ysize; j++ ) { for ( i = 0; i < xsize; i++ ) { if ( rgb_max < *r_index ) { rgb_max = *r_index; } r_index = r_index + 1; if ( rgb_max < *g_index ) { rgb_max = *g_index; } g_index = g_index + 1; if ( rgb_max < *b_index ) { rgb_max = *b_index; } b_index = b_index + 1; } } /* Write the header. */ error = ppma_write_header ( file_out, file_out_name, xsize, ysize, rgb_max ); if ( error ) { printf ( "\n" ); printf ( "PPMA_WRITE - Fatal error!\n" ); printf ( " PPMA_WRITE_HEADER failed.\n" ); return 1; } /* Write the data. */ error = ppma_write_data ( file_out, xsize, ysize, r, g, b ); if ( error ) { printf ( "\n" ); printf ( "PPMA_WRITE - Fatal error!\n" ); printf ( " PPMA_WRITE_DATA failed.\n" ); return 1; } /* Close the file. */ fclose ( file_out ); return 0; } /******************************************************************************/ int ppma_write_data ( FILE *file_out, int xsize, int ysize, int *r, int *g, int *b ) /******************************************************************************/ /* Purpose: PPMA_WRITE_DATA writes the data for an ASCII portable pixel map file. Licensing: This code is distributed under the GNU LGPL license. Modified: 28 February 2003 Author: John Burkardt Parameters: Input, ofstream &FILE_OUT, a pointer to the file to contain the ASCII portable pixel map data. Input, int XSIZE, YSIZE, the number of rows and columns of data. Input, int *R, *G, *B, the arrays of XSIZE by YSIZE data values. Output, int PPMA_WRITE_DATA, is true, if an error was detected, or false, if the data was written. */ { int *b_index; int *g_index; int i; int j; int *r_index; int rgb_num; r_index = r; g_index = g; b_index = b; rgb_num = 0; for ( j = 0; j < ysize; j++ ) { for ( i = 0; i < xsize; i++ ) { fprintf ( file_out, "%d %d %d", *r_index, *g_index, *b_index ); rgb_num = rgb_num + 3; r_index = r_index + 1; g_index = g_index + 1; b_index = b_index + 1; if ( rgb_num % 12 == 0 || i == xsize - 1 || rgb_num == 3 * xsize * ysize ) { fprintf ( file_out, "\n" ); } else { fprintf ( file_out, " " ); } } } return 0; } /******************************************************************************/ int ppma_write_header ( FILE *file_out, char *file_out_name, int xsize, int ysize, int rgb_max ) /******************************************************************************/ /* Purpose: PPMA_WRITE_HEADER writes the header of an ASCII portable pixel map file. Licensing: This code is distributed under the GNU LGPL license. Modified: 28 February 2003 Author: John Burkardt Parameters: Input, FILE *FILE_OUT, a pointer to the file to contain the ASCII portable pixel map data. Input, char *FILE_OUT_NAME, the name of the file. Input, int XSIZE, YSIZE, the number of rows and columns of data. Input, int RGB_MAX, the maximum RGB value. Output, int PPMA_WRITE_HEADER, is true, if an error was detected, or false, if the header was written. */ { fprintf ( file_out, "P3\n" ); fprintf ( file_out, "# %s created by ppma_write.c.\n", file_out_name ); fprintf ( file_out, "%d %d\n", xsize, ysize ); fprintf ( file_out, "%d\n", rgb_max ); return 0; }