Index: sm/trunk/src/c/modules/HoleFillerx/HoleFillerCrippen.c =================================================================== --- /issm/trunk/src/c/modules/HoleFillerx/HoleFillerCrippen.c (revision 9312) +++ (revision ) @@ -1,1296 +1,0 @@ -/* Hole Filler Version 2003.04.22 */ -/* Grow hole edges one pixel at a time by interpolation. */ -/* */ -/* Only fills void pixels that are next to non-void pixels. */ -/* Interpolates in eight directions using weighted average. */ -/* Iterates on the entire image until all voids are filled. */ -/* Thus grows edges by interpolation - best of both worlds. */ -/* Smoothing option is slow but effective (two-byte data only). */ - -/* There is a lot of experimental code here. */ -/* The program was intended to do two iterations in each pass */ -/* through the code and the screen output indicates that. */ -/* Input A --> Output B becomes Input B --> Output A in memory. */ -/* For experimental purposes this was bypassed so that the code */ -/* could be written just once. */ -/* Instead, for now, Output B is copied to Input A for the next */ -/* iteration. */ - - -#include -#include -#include - - -#include -#include - -int number_input(); -int yes_no_input(); - -main(void) -{ - FILE *fp1, *fp2; - char file1name[400], file2name[400], pause[400]; - unsigned long filesize; - long lines, samps, infinit, temp; - long i, j, k, ii, jj, bytes; - long iii, jjj; - long count, test, prevlines, prevpix, test1, test2, thresh; - long bell, cycles, leftovers, smooth; - float howlong; - float nsteps, ssteps, wsteps, esteps; - float nwsteps, nesteps, swsteps, sesteps; - float north, south, east, west; - float nw, sw, ne, se, elev, range; - float nnesteps, enesteps, esesteps, ssesteps; - float sswsteps, wswsteps, wnwsteps, nnwsteps; - float nne, ene, ese, sse; - float ssw, wsw, wnw, nnw; - float sum; - - time_t t1, t2; - - unsigned char **dns1 , **dns2 , **dns3 , *buffy ; - unsigned short **dns01, **dns02, **dns03, *buffy2; - - printf( "\n" ); - - printf( ">>> Hole Filler Version 2003.04.22 <<<\n" ); - printf( ">>> Grow hole edges one pixel at a time by interpolation. <<<\n" ); - printf( ">>> <<<\n" ); - printf( ">>> Only fills void pixels that are next to non-void pixels. <<<\n" ); - printf( ">>> Interpolates in eight directions using weighted average. <<<\n" ); - printf( ">>> Iterates on the entire image until all voids are filled. <<<\n" ); - printf( ">>> Thus grows edges by interpolation - best of both worlds. <<<\n" ); - printf( ">>> Smoothing option is slow but effective (two-byte data only). <<<\n" ); - - retry1: - printf( "\nEnter INPUT IMAGE file NAME: " ); - gets(file1name); - if ( (fp1 = fopen( file1name, "rb")) == NULL ) - { printf( "*** No such file found! *** \a\n" ); - goto retry1; } - - fseek( fp1, 0, SEEK_END ); - filesize = ftell( fp1 ); - fseek( fp1, 0, SEEK_SET ); - - count = 0; - - printf( "\nData Type ( Default = ONE-BYTE ), type '2' for two-byte data: "); - bytes = number_input(); if ( bytes > 1 ) bytes = 2; else bytes = 1; - if ( bytes == 1 ) printf( "ONE-BYTE data specified.\n"); - else { printf( "TWO-BYTE data specified.\n"); goto twobyte; } - - - retry2: - do { printf( "\nNumber of SAMPS in image: "); - samps = number_input(); - if (samps < 1) printf( "*** Must be positive integer. *** \a\n" ); - } while ( samps < 1 ); - - if ( filesize % samps ) - { - printf( "*** Number of SAMPS is INCORRECT for file. *** \n" ); - printf( "*** It has %ld pixels. ", filesize ); - printf( "That is %6.5f lines. ***\n\a", (float)filesize / samps ); - goto retry2; - } - - lines = filesize / samps; - printf( "Number of LINES in image = %ld \n", lines ); - infinit = lines + samps; - - /* No smoothing yet for one-byte. - - retry4: - printf( "\nSMOOTH the output with a boxfilter? ( Default = Yes, else type 'no'): "); - smooth = yes_no_input(); - if ( smooth == 0 ) printf("Data patches will be SMOOTHED. \n"); - if ( smooth == 1 ) printf("Data patches will NOT be smoothed. \n"); - - */ - - retry3: - printf( "\nEnter OUTPUT IMAGE file NAME: " ); - gets(file2name); - if ( strcmp(file1name,file2name) == 0 ) - { - printf("*** Cannot overwrite input file! ***\n\a"); - goto retry3; - } - if ( strlen( file2name ) == 0 ) { printf("\a"); goto retry3; } - if ( (fp2 = fopen( file2name, "rb")) != NULL ) - { - fclose( fp2 ); - printf( "*** File already EXISTS! *** \a\n\n" ); - printf( "Hit 'return' to REPLACE, type 'no' to rename: "); - if ( yes_no_input() > 0 ) goto retry3; - } - howlong = filesize / 44347.0 / 60 ; /* 44347 = empirical speed factor */ - - - - - - printf("\nPull down menu to 'Quit', else hit 'return':"); gets(pause); - - time(&t1); printf( "\nStart "); printf( ctime(&t1) ); - - fp2 = fopen( file2name, "wb" ); - - - dns1 = (unsigned char **) malloc( lines * sizeof( char * ) ); - if (dns1 == NULL) { printf("Not enough RAM.\a\n"); goto end; } - - for ( i = 0; i < lines; i ++ ) /* Get input RAM */ - { - dns1[ i ] = ( unsigned char * ) malloc( samps ); - if (dns1[ i ] == NULL) { printf("Not enough RAM.\a\n"); goto end; } - } - - - /*vvvvvvvvvvvv Start Buffering of Input vvvvvvvvvvvv*/ - - buffy = (unsigned char *) malloc( samps*100 ); - if (buffy == NULL) { printf("Not enough RAM.\a\n"); goto end; } - - cycles = lines / 100; - for ( k = 0; k < cycles; k++ ) - { - fread( buffy, 1, 100*samps, fp1 ); - for ( i = 0; i < 100; i++ ) - { - prevlines = ( k * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - dns1 [ prevlines ][ j ] = buffy [ prevpix + j ]; - } - } - - leftovers = lines - ( cycles * 100 ); - fread( buffy, 1, samps * leftovers, fp1 ); - for ( i = 0; i < leftovers; i++ ) - { - prevlines = ( cycles * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - dns1 [ prevlines ][ j ] = buffy [ prevpix + j ]; - } - - /*^^^^^^^^^^^^^ End Buffering of Input ^^^^^^^^^^^^^*/ - - - dns2 = (unsigned char **) malloc( lines * sizeof( char * ) ); - if (dns2 == NULL) { printf("Not enough RAM.\a\n"); goto end; } - - for ( i = 0; i < lines; i ++ ) /* Get output RAM */ - { - dns2[ i ] = ( unsigned char * ) malloc( samps ); - if (dns2[ i ] == NULL) { printf("Not enough RAM.\a\n"); goto end; } - } - - for ( i = 0; i < lines; i++ ) - for ( j = 0; j < samps; j++ ) - dns2[i][j] = 0; - - -/***************** Cycling of Run Pairs *********************/ - - count += 2; - printf( "\n" ); - printf( "Iterations %5ld and %5ld.", count-1, count ); - fflush( stdout ); - goto afterfirst1; - - - again1: - - count += 2; - - printf( "\b\b\b\b\b\b\b\b\b\b\b\b\b\b" ); - printf( "\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b" ); - printf( "Iterations %5ld and %5ld.", count-1, count ); - fflush( stdout ); - - afterfirst1: - - -/***************** FIRST RUN *********************/ - -/* - fprint ( "First Application: " ); time(&t2); printf( ctime(&t2) ); -*/ - - for ( i = 0; i < lines; i++ ) - { - for ( j = 0; j < samps; j++ ) - { - if ( dns1[i][j] > 0 ) { dns2[i][j] = dns1[i][j]; continue; } /* If not a void pixel */ - - test = 0; - for ( ii = i-1; ii < i+2; ii++ ) /* Find it to be a void edge */ - { - for ( jj = j-1; jj < j+2; jj++ ) - { - if (( ii < 0 ) || ( ii >= lines )) continue; - if (( jj < 0 ) || ( jj >= samps )) continue; - if ( dns1[ii][jj] > 0 ) test = 1; - } - } - if ( test == 0 ) { dns2[i][j] = 0; continue; } /* If not a void edge */ - - - /* For void edge pixels: */ - - nsteps = 0.0; ssteps = 0.0; esteps = 0.0; wsteps = 0.0; - nwsteps = 0.0; nesteps = 0.0; swsteps = 0.0; sesteps = 0.0; - north = 0.0; south = 0.0; east = 0.0; west = 0.0; - nw = 0.0; ne = 0.0; sw = 0.0; se = 0.0; - - - - for ( ii = i - 1; ii >= 0; ii-- ) /* North */ - { - if ( ii <= 0 ) { north = 0; nsteps = 0; break; } - if ( dns1[ii][j] > 0 ) { north = dns1[ii][j]; nsteps = i - ii; break; } - } - - for ( ii = i + 1; ii < lines; ii++ ) /* South */ - { - if ( ii >= lines-1 ) { south = 0; ssteps = 0; break; } - if ( dns1[ii][j] > 0 ) { south = dns1[ii][j]; ssteps = ii - i; break; } - } - - for ( jj = j - 1; jj >= 0; jj-- ) /* West */ - { - if ( jj <= 0 ) { west = 0; wsteps = 0; break; } - if ( dns1[i][jj] > 0 ) { west = dns1[i][jj]; wsteps = j - jj; break; } - } - - for ( jj = j + 1; jj < samps; jj++ ) /* East */ - { - if ( jj >= samps-1 ) { east = 0; esteps = 0; break; } - if ( dns1[i][jj] > 0 ) { east = dns1[i][jj]; esteps = jj - j; break; } - } - - - - /* Southeast */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j + k; - if (( ii == lines ) || ( jj == samps )) { se = 0; sesteps = 0; break; } - if ( dns1[ii][jj] > 0 ) { se = dns1[ii][jj]; sesteps = 1.4142 * k; break; } - } - - - /* Northeast */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j + k; - if (( ii == -1 ) || ( jj == samps )) { ne = 0; nesteps = 0; break; } - if ( dns1[ii][jj] > 0 ) { ne = dns1[ii][jj]; nesteps = 1.4142 * k; break; } - } - - /* Northwest */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j - k; - if (( ii == -1 ) || ( jj == -1 )) { nw = 0; nwsteps = 0; break; } - if ( dns1[ii][jj] > 0 ) { nw = dns1[ii][jj]; nwsteps = 1.4142 * k; break; } - } - - - /* Southwest */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j - k; - if (( ii == lines ) || ( jj == -1 )) { sw = 0; swsteps = 0; break; } - if ( dns1[ii][jj] > 0 ) { sw = dns1[ii][jj]; swsteps = 1.4142 * k; break; } - } - - - elev = 0; range = 0; - if ( nsteps > 0.5 ) { elev += north / nsteps; range += 1.00 / nsteps; } - if ( ssteps > 0.5 ) { elev += south / ssteps; range += 1.00 / ssteps; } - if ( esteps > 0.5 ) { elev += east / esteps; range += 1.00 / esteps; } - if ( wsteps > 0.5 ) { elev += west / wsteps; range += 1.00 / wsteps; } - if ( nwsteps > 0.5 ) { elev += nw / nwsteps; range += 1.00 / nwsteps; } - if ( nesteps > 0.5 ) { elev += ne / nesteps; range += 1.00 / nesteps; } - if ( swsteps > 0.5 ) { elev += sw / swsteps; range += 1.00 / swsteps; } - if ( sesteps > 0.5 ) { elev += se / sesteps; range += 1.00 / sesteps; } - - - temp = ( elev / range ) + 0.5 ; - if ( temp > 255 ) temp = 255; - if ( temp < 0 ) temp = 0; - - dns2[i][j] = temp ; - - } - } - - -/***************** SECOND RUN *********************/ - -/* -fprintf ( stdout, "Second Application: " ); time(&t2); printf( ctime(&t2) ); -*/ - - for ( i = 0; i < lines; i++ ) - { - for ( j = 0; j < samps; j++ ) - { - if ( dns2[i][j] > 0 ) { dns1[i][j] = dns2[i][j]; continue; } /* If not a void pixel */ - - test = 0; - for ( ii = i-1; ii < i+2; ii++ ) /* Find it to be a void edge */ - { - for ( jj = j-1; jj < j+2; jj++ ) - { - if (( ii < 0 ) || ( ii >= lines )) continue; - if (( jj < 0 ) || ( jj >= samps )) continue; - if ( dns2[ii][jj] > 0 ) test = 1; - } - } - if ( test == 0 ) { dns1[i][j] = 0; continue; } /* If not a void edge */ - - - /* For void edge pixels: */ - - nsteps = 0.0; ssteps = 0.0; esteps = 0.0; wsteps = 0.0; - nwsteps = 0.0; nesteps = 0.0; swsteps = 0.0; sesteps = 0.0; - north = 0.0; south = 0.0; east = 0.0; west = 0.0; - nw = 0.0; ne = 0.0; sw = 0.0; se = 0.0; - - - - for ( ii = i - 1; ii >= 0; ii-- ) /* North */ - { - if ( ii <= 0 ) { north = 0; nsteps = 0; break; } - if ( dns2[ii][j] > 0 ) { north = dns2[ii][j]; nsteps = i - ii; break; } - } - - for ( ii = i + 1; ii < lines; ii++ ) /* South */ - { - if ( ii >= lines-1 ) { south = 0; ssteps = 0; break; } - if ( dns2[ii][j] > 0 ) { south = dns2[ii][j]; ssteps = ii - i; break; } - } - - for ( jj = j - 1; jj >= 0; jj-- ) /* West */ - { - if ( jj <= 0 ) { west = 0; wsteps = 0; break; } - if ( dns2[i][jj] > 0 ) { west = dns2[i][jj]; wsteps = j - jj; break; } - } - - for ( jj = j + 1; jj < samps; jj++ ) /* East */ - { - if ( jj >= samps-1 ) { east = 0; esteps = 0; break; } - if ( dns2[i][jj] > 0 ) { east = dns2[i][jj]; esteps = jj - j; break; } - } - - - - /* Southeast */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j + k; - if (( ii == lines ) || ( jj == samps )) { se = 0; sesteps = 0; break; } - if ( dns2[ii][jj] > 0 ) { se = dns2[ii][jj]; sesteps = 1.4142 * k; break; } - } - - - /* Northeast */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j + k; - if (( ii == -1 ) || ( jj == samps )) { ne = 0; nesteps = 0; break; } - if ( dns2[ii][jj] > 0 ) { ne = dns2[ii][jj]; nesteps = 1.4142 * k; break; } - } - - /* Northwest */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j - k; - if (( ii == -1 ) || ( jj == -1 )) { nw = 0; nwsteps = 0; break; } - if ( dns2[ii][jj] > 0 ) { nw = dns2[ii][jj]; nwsteps = 1.4142 * k; break; } - } - - - /* Southwest */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j - k; - if (( ii == lines ) || ( jj == -1 )) { sw = 0; swsteps = 0; break; } - if ( dns2[ii][jj] > 0 ) { sw = dns2[ii][jj]; swsteps = 1.4142 * k; break; } - } - - - elev = 0; range = 0; - if ( nsteps > 0.5 ) { elev += north / nsteps; range += 1.00 / nsteps; } - if ( ssteps > 0.5 ) { elev += south / ssteps; range += 1.00 / ssteps; } - if ( esteps > 0.5 ) { elev += east / esteps; range += 1.00 / esteps; } - if ( wsteps > 0.5 ) { elev += west / wsteps; range += 1.00 / wsteps; } - if ( nwsteps > 0.5 ) { elev += nw / nwsteps; range += 1.00 / nwsteps; } - if ( nesteps > 0.5 ) { elev += ne / nesteps; range += 1.00 / nesteps; } - if ( swsteps > 0.5 ) { elev += sw / swsteps; range += 1.00 / swsteps; } - if ( sesteps > 0.5 ) { elev += se / sesteps; range += 1.00 / sesteps; } - - - temp = ( elev / range ) + 0.5 ; - if ( temp > 255 ) temp = 255; - if ( temp < 0 ) temp = 0; - - dns1[i][j] = temp ; - - } - } - - -/***************** CHECK FOR REMAINING ZEROES *********************/ - - - for ( i = 0; i < lines; i++ ) - for ( j = 0; j < samps; j++ ) - if ( dns1 [i][j] == 0 ) goto again1; - - - - /*vvvvvvvvvvvv Start Buffering of Output vvvvvvvvvvvv*/ - - cycles = lines / 100; - for ( k = 0; k < cycles; k++ ) - { - for ( i = 0; i < 100; i++ ) - { - prevlines = ( k * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - buffy [ prevpix + j ] = dns1 [ prevlines ][ j ]; - } - fwrite( buffy, 1, 100*samps, fp2 ); - } - - leftovers = lines - ( cycles * 100 ); - for ( i = 0; i < leftovers; i++ ) - { - prevlines = ( cycles * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - buffy [ prevpix + j ] = dns1 [ prevlines ][ j ]; - } - fwrite( buffy, 1, samps * leftovers, fp2 ); - - - /*^^^^^^^^^^^^^ End Buffering of Output ^^^^^^^^^^^^^*/ - - - fclose( fp1 ); fclose( fp2 ); - - time(&t2); printf ( "\n\nEnd "); printf( ctime(&t2) ); - printf( "\n Total time: "); - printf( "%4.2f minutes. \n", difftime(t2,t1)/60 ); - goto end; - - - - - - - -/*============================= TWO-BYTE VERSION ================================*/ - - twobyte: - - retry02: - do { printf( "\nNumber of SAMPS in image: "); - samps = number_input(); - if (samps < 1) printf( "*** Must be positive integer. *** \a\n" ); - } while ( samps < 1 ); - - if ( filesize % samps ) - { - printf( "*** Number of SAMPS is INCORRECT for file. *** \n" ); - printf( "*** It has %ld pixels. ", filesize ); - printf( "That is %6.5f lines. ***\n\a", (float)filesize / samps ); - goto retry02; - } - - lines = filesize / ( samps * sizeof(short) ); - printf( "Number of LINES in image = %ld \n", lines ); - infinit = lines + samps; - - - - retry04: - printf( "\nSMOOTH the output with a boxfilter? ( Default = Yes, else type 'no'): "); - smooth = yes_no_input(); - if ( smooth == 0 ) printf("Data patches will be SMOOTHED. \n"); - if ( smooth == 1 ) printf("Data patches will NOT be smoothed. \n"); - - - - - retry03: - printf( "\nEnter OUTPUT IMAGE file NAME: " ); - gets(file2name); - if ( strcmp(file1name,file2name) == 0 ) - { - printf("*** Cannot overwrite input file! ***\n\a"); - goto retry03; - } - if ( strlen( file2name ) == 0 ) { printf("\a"); goto retry03; } - if ( (fp2 = fopen( file2name, "rb")) != NULL ) - { - fclose( fp2 ); - printf( "*** File already EXISTS! *** \a\n\n" ); - printf( "Hit 'return' to REPLACE, type 'no' to rename: "); - if ( yes_no_input() > 0 ) goto retry03; - } - howlong = filesize / 44347.0 / 60 ; /* 44347 = empirical speed factor */ - - - - - printf("\nPull down menu to 'Quit', else hit 'return':"); gets(pause); - - time(&t1); printf( "\nStart "); printf( ctime(&t1) ); - - fp2 = fopen( file2name, "wb" ); - - - dns01 = (unsigned short **) malloc( lines * sizeof( short * ) ); - if (dns01 == NULL) { printf("Not enough RAM.\a\n"); goto end; } - - for ( i = 0; i < lines; i ++ ) /* Get input file */ - { - dns01[ i ] = ( unsigned short * ) malloc( samps * sizeof(short) ); - if (dns01[ i ] == NULL) { printf("Not enough RAM.\a\n"); goto end; } - } - - - /*vvvvvvvvvvvv Start Buffering of Input vvvvvvvvvvvv*/ - - buffy2 = (unsigned short *) malloc( samps*100 * sizeof(short) ); - if (buffy2 == NULL) { printf("Not enough RAM.\a\n"); goto end; } - - cycles = lines / 100; - for ( k = 0; k < cycles; k++ ) - { - fread( buffy2, sizeof(short), 100*samps, fp1 ); - for ( i = 0; i < 100; i++ ) - { - prevlines = ( k * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - dns01 [ prevlines ][ j ] = buffy2 [ prevpix + j ]; - } - } - - leftovers = lines - ( cycles * 100 ); - fread( buffy2, sizeof(short), samps * leftovers, fp1 ); - for ( i = 0; i < leftovers; i++ ) - { - prevlines = ( cycles * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - dns01 [ prevlines ][ j ] = buffy2 [ prevpix + j ]; - } - - /*^^^^^^^^^^^^^ End Buffering of Input ^^^^^^^^^^^^^*/ - - - dns02 = (unsigned short **) malloc( lines * sizeof( short * ) ); - if (dns02 == NULL) { printf("Not enough RAM.\a\n"); goto end; } - - for ( i = 0; i < lines; i ++ ) /* Get output RAM */ - { - dns02[ i ] = ( unsigned short * ) malloc( samps * sizeof(short) ); - if (dns02[ i ] == NULL) { printf("Not enough RAM.\a\n"); goto end; } - } - - -/********** NEW: Eliminate loner pixels ***************/ - -/* - for ( i = 0; i < lines; i++ ) - { - for ( j = 0; j < samps; j++ ) - { - dns02[i][j] = dns01[i][j]; - test = 0; - if ((( i == 0 ) || ( dns01[i-1][j] < 1 ) ) && ( ( i == lines-1 ) || ( dns01[i+1][j] < 1 ))) test += 1; - if ((( j == 0 ) || ( dns01[i][j-1] < 1 ) ) && ( ( j == samps-1 ) || ( dns01[i][j+1] < 1 ))) test += 1; - if ( test > 0 ) dns02[i][j] = 0; - } - } - - for ( i = 0; i < lines; i++ ) - { - for ( j = 0; j < samps; j++ ) - { - dns01[i][j] = dns02[i][j]; - } - } -*/ - - -/***************** Cycling of Run Pairs *********************/ - - - - count += 2; - printf( "\n" ); - printf( "Iterations %5ld and %5ld.", count-1, count ); - fflush( stdout ); - goto afterfirst2; - - - again2: - - count += 2; - - printf( "\b\b\b\b\b\b\b\b\b\b\b\b\b\b" ); - printf( "\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b" ); - printf( "Iterations %5ld and %5ld.", count-1, count ); - fflush( stdout ); - - afterfirst2: - - - - -/***************** FIRST RUN *********************/ - -/* -fprintf ( stdout, "First Application: " ); time(&t2); printf( ctime(&t2) ); -*/ - - for ( i = 0; i < lines; i++ ) - { - for ( j = 0; j < samps; j++ ) - { - if ( dns01[i][j] > 0 ) { dns02[i][j] = dns01[i][j]; continue; } /* If not a void pixel */ - - - /* Find it to be a void edge */ - -/* Sides Only Test1 */ -/* - test1 = 0; - if (( i > 0 ) && ( dns01[i-1][j] > 0 )) test1 += 1; - if (( i < lines-1 ) && ( dns01[i+1][j] > 0 )) test1 += 1; - if (( j > 0 ) && ( dns01[i][j-1] > 0 )) test1 += 1; - if (( j < samps-1 ) && ( dns01[i][j+1] > 0 )) test1 += 1; - - if ( test1 < 1 ) { dns02[i][j] = 0; continue; } -*/ -/* Diagnals Only Test2 */ - -/* test2 = 0; - if (( i > 0 ) && ( j > 0 ) && ( dns01[i-1][j-1] > 0 )) test2 += 1; - if (( i < lines-1 ) && ( j < samps-1 ) && ( dns01[i+1][j+1] > 0 )) test2 += 1; - if (( i > 0 ) && ( j < samps-1 ) && ( dns01[i-1][j+1] > 0 )) test2 += 1; - if (( i < lines-1 ) && ( j > 0 ) && ( dns01[i+1][j-1] > 0 )) test2 += 1; - - if ( test2 < 2 ) { dns02[i][j] = 0; continue; } - - if (( test1 < 1 ) && ( test2 < 2 )) { dns02[i][j] = 0; continue; } -*/ - - - - - - -/* Count All Neighbors Test */ - - test = 0; - for ( ii = i-1; ii < i+2; ii++ ) - { - for ( jj = j-1; jj < j+2; jj++ ) - { - if (( ii < 0 ) || ( ii >= lines )) continue; - if (( jj < 0 ) || ( jj >= samps )) continue; - if ( dns01[ii][jj] > 0 ) test += 1; - } - } - if ( test < 3 ) { dns02[i][j] = 0; continue; } - - - - - - - - - - /* For void edge pixels: */ - - nsteps = 0.0; ssteps = 0.0; esteps = 0.0; wsteps = 0.0; - nwsteps = 0.0; nesteps = 0.0; swsteps = 0.0; sesteps = 0.0; - north = 0.0; south = 0.0; east = 0.0; west = 0.0; - nw = 0.0; ne = 0.0; sw = 0.0; se = 0.0; - - nnesteps = 0.0; enesteps = 0.0; esesteps = 0.0; ssesteps = 0.0; - sswsteps = 0.0; wswsteps = 0.0; wnwsteps = 0.0; nnwsteps = 0.0; - nne = 0.0; ene = 0.0; ese = 0.0; sse = 0.0; - ssw = 0.0; wsw = 0.0; wnw = 0.0; nnw = 0.0; - - -/** NSEW **/ - - for ( ii = i - 1; ii >= 0; ii-- ) /* North */ - { - if ( ii <= 0 ) { north = 0; nsteps = 0; break; } - if ( dns01[ii][j] > 0 ) { north = dns01[ii][j]; nsteps = i - ii; break; } - } - - for ( ii = i + 1; ii < lines; ii++ ) /* South */ - { - if ( ii >= lines-1 ) { south = 0; ssteps = 0; break; } - if ( dns01[ii][j] > 0 ) { south = dns01[ii][j]; ssteps = ii - i; break; } - } - - for ( jj = j - 1; jj >= 0; jj-- ) /* West */ - { - if ( jj <= 0 ) { west = 0; wsteps = 0; break; } - if ( dns01[i][jj] > 0 ) { west = dns01[i][jj]; wsteps = j - jj; break; } - } - - for ( jj = j + 1; jj < samps; jj++ ) /* East */ - { - if ( jj >= samps-1 ) { east = 0; esteps = 0; break; } - if ( dns01[i][jj] > 0 ) { east = dns01[i][jj]; esteps = jj - j; break; } - } - - -/** Diagonals **/ - - /* Southeast */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j + k; - if (( ii == lines ) || ( jj == samps )) { se = 0; sesteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { se = dns01[ii][jj]; sesteps = 1.4142 * k; break; } - } - - - /* Northeast */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j + k; - if (( ii == -1 ) || ( jj == samps )) { ne = 0; nesteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { ne = dns01[ii][jj]; nesteps = 1.4142 * k; break; } - } - - /* Northwest */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j - k; - if (( ii == -1 ) || ( jj == -1 )) { nw = 0; nwsteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { nw = dns01[ii][jj]; nwsteps = 1.4142 * k; break; } - } - - - /* Southwest */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j - k; - if (( ii == lines ) || ( jj == -1 )) { sw = 0; swsteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { sw = dns01[ii][jj]; swsteps = 1.4142 * k; break; } - } - - -/** Other 8 **/ - - /* NNE */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - (k*2); jj = j + k; - if (( ii <= -1 ) || ( jj >= samps )) { nne = 0; nnesteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { nne = dns01[ii][jj]; nnesteps = 2.2361 * k; break; } - } - - - /* ENE */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j + (k*2); - if (( ii <= -1 ) || ( jj >= samps )) { ene = 0; enesteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { ene = dns01[ii][jj]; enesteps = 2.2361 * k; break; } - } - - /* ESE */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j + (k*2); - if (( ii >= lines ) || ( jj >= samps )) { ese = 0; esesteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { ese = dns01[ii][jj]; esesteps = 2.2361 * k; break; } - } - - - /* SSE */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + (k*2); jj = j + k; - if (( ii >= lines ) || ( jj >= samps )) { sse = 0; ssesteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { sse = dns01[ii][jj]; ssesteps = 2.2361 * k; break; } - } - - - - /* SSW */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + (k*2); jj = j - k; - if (( ii >= lines ) || ( jj <= -1 )) { ssw = 0; sswsteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { ssw = dns01[ii][jj]; sswsteps = 2.2361 * k; break; } - } - - - /* WSW */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j - (k*2); - if (( ii >= lines ) || ( jj <= -1 )) { wsw = 0; wswsteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { wsw = dns01[ii][jj]; wswsteps = 2.2361 * k; break; } - } - - /* WNW */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j - (k*2); - if (( ii <= -1 ) || ( jj <= -1 )) { wnw = 0; wnwsteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { wnw = dns01[ii][jj]; wnwsteps = 2.2361 * k; break; } - } - - - /* NNW */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - (k*2); jj = j - k; - if (( ii <= -1 ) || ( jj <= -1 )) { nnw = 0; nnwsteps = 0; break; } - if ( dns01[ii][jj] > 0 ) { nnw = dns01[ii][jj]; nnwsteps = 2.2361 * k; break; } - } - - - -/* Avoid fills that use too many distant points */ -/* -thresh = 0; - if ( nsteps > 5 ) thresh += 1; - if ( ssteps > 5 ) thresh += 1; - if ( esteps > 5 ) thresh += 1; - if ( wsteps > 5 ) thresh += 1; - - if ( nwsteps > 5 ) thresh += 1; - if ( nesteps > 5 ) thresh += 1; - if ( swsteps > 5 ) thresh += 1; - if ( sesteps > 5 ) thresh += 1; - - if ( nnesteps > 5 ) thresh += 1; - if ( enesteps > 5 ) thresh += 1; - if ( esesteps > 5 ) thresh += 1; - if ( ssesteps > 5 ) thresh += 1; - if ( sswsteps > 5 ) thresh += 1; - if ( wswsteps > 5 ) thresh += 1; - if ( wnwsteps > 5 ) thresh += 1; - if ( nnwsteps > 5 ) thresh += 1; -if ( thresh > 14 ) { dns02[i][j] = 0; continue; } -*/ - - - -/* TRY SQUARES -- BAD! */ -/* - nsteps = pow(nsteps, 0.5); - ssteps = pow(ssteps, 0.5); - esteps = pow(esteps, 0.5); - wsteps = pow(wsteps, 0.5); - nwsteps = pow(nwsteps, 0.5); - nesteps = pow(nesteps, 0.5); - swsteps = pow(swsteps, 0.5); - sesteps = pow(sesteps, 0.5); - nnesteps = pow(nnesteps, 0.5); - enesteps = pow(enesteps, 0.5); - esesteps = pow(esesteps, 0.5); - ssesteps = pow(ssesteps, 0.5); - sswsteps = pow(sswsteps, 0.5); - wswsteps = pow(wswsteps, 0.5); - wnwsteps = pow(wnwsteps, 0.5); - nnwsteps = pow(nnwsteps, 0.5); -*/ - - - - elev = 0; range = 0; -/*NSEW*/ - if ( nsteps > 0.5 ) { elev += north / nsteps; range += 1.00 / nsteps; } - if ( ssteps > 0.5 ) { elev += south / ssteps; range += 1.00 / ssteps; } - if ( esteps > 0.5 ) { elev += east / esteps; range += 1.00 / esteps; } - if ( wsteps > 0.5 ) { elev += west / wsteps; range += 1.00 / wsteps; } -/*Diag*/ - if ( nwsteps > 0.5 ) { elev += nw / nwsteps; range += 1.00 / nwsteps; } - if ( nesteps > 0.5 ) { elev += ne / nesteps; range += 1.00 / nesteps; } - if ( swsteps > 0.5 ) { elev += sw / swsteps; range += 1.00 / swsteps; } - if ( sesteps > 0.5 ) { elev += se / sesteps; range += 1.00 / sesteps; } - -/*Other 8*/ - if ( nnesteps > 0.5 ) { elev += nne / nnesteps; range += 1.00 / nnesteps; } - if ( enesteps > 0.5 ) { elev += ene / enesteps; range += 1.00 / enesteps; } - if ( esesteps > 0.5 ) { elev += ese / esesteps; range += 1.00 / esesteps; } - if ( ssesteps > 0.5 ) { elev += sse / ssesteps; range += 1.00 / ssesteps; } - if ( sswsteps > 0.5 ) { elev += ssw / sswsteps; range += 1.00 / sswsteps; } - if ( wswsteps > 0.5 ) { elev += wsw / wswsteps; range += 1.00 / wswsteps; } - if ( wnwsteps > 0.5 ) { elev += wnw / wnwsteps; range += 1.00 / wnwsteps; } - if ( nnwsteps > 0.5 ) { elev += nnw / nnwsteps; range += 1.00 / nnwsteps; } - - - temp = ( elev / range ) + 0.5 ; - if ( temp > 10000 ) temp = 10000; - if ( temp < 0 ) temp = 0; - - dns02[i][j] = temp ; - - } - } - - - - -goto there; - - - - - -/***************** SECOND RUN *********************/ - -fprintf ( stdout, "Second Application: " ); time(&t2); printf( ctime(&t2) ); - - - for ( i = 0; i < lines; i++ ) - { - for ( j = 0; j < samps; j++ ) - { - if ( dns02[i][j] > 0 ) { dns01[i][j] = dns02[i][j]; continue; } /* If not a void pixel */ - - test = 0; - for ( ii = i-1; ii < i+2; ii++ ) /* Find it to be a void edge */ - { - for ( jj = j-1; jj < j+2; jj++ ) - { - if (( ii < 0 ) || ( ii >= lines )) continue; - if (( jj < 0 ) || ( jj >= samps )) continue; - if ( dns02[ii][jj] > 0 ) test = 1; - } - } - if ( test == 0 ) { dns01[i][j] = 0; continue; } /* If not a void edge */ - - - /* For void edge pixels: */ - - nsteps = 0.0; ssteps = 0.0; esteps = 0.0; wsteps = 0.0; - nwsteps = 0.0; nesteps = 0.0; swsteps = 0.0; sesteps = 0.0; - north = 0.0; south = 0.0; east = 0.0; west = 0.0; - nw = 0.0; ne = 0.0; sw = 0.0; se = 0.0; - - - - for ( ii = i - 1; ii >= 0; ii-- ) /* North */ - { - if ( ii <= 0 ) { north = 0; nsteps = 0; break; } - if ( dns02[ii][j] > 0 ) { north = dns02[ii][j]; nsteps = i - ii; break; } - } - - for ( ii = i + 1; ii < lines; ii++ ) /* South */ - { - if ( ii >= lines-1 ) { south = 0; ssteps = 0; break; } - if ( dns02[ii][j] > 0 ) { south = dns02[ii][j]; ssteps = ii - i; break; } - } - - for ( jj = j - 1; jj >= 0; jj-- ) /* West */ - { - if ( jj <= 0 ) { west = 0; wsteps = 0; break; } - if ( dns02[i][jj] > 0 ) { west = dns02[i][jj]; wsteps = j - jj; break; } - } - - for ( jj = j + 1; jj < samps; jj++ ) /* East */ - { - if ( jj >= samps-1 ) { east = 0; esteps = 0; break; } - if ( dns02[i][jj] > 0 ) { east = dns02[i][jj]; esteps = jj - j; break; } - } - - - - /* Southeast */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j + k; - if (( ii == lines ) || ( jj == samps )) { se = 0; sesteps = 0; break; } - if ( dns02[ii][jj] > 0 ) { se = dns02[ii][jj]; sesteps = 1.4142 * k; break; } - } - - - /* Northeast */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j + k; - if (( ii == -1 ) || ( jj == samps )) { ne = 0; nesteps = 0; break; } - if ( dns02[ii][jj] > 0 ) { ne = dns02[ii][jj]; nesteps = 1.4142 * k; break; } - } - - /* Northwest */ - for ( k = 1; k < infinit; k++ ) - { - ii = i - k; jj = j - k; - if (( ii == -1 ) || ( jj == -1 )) { nw = 0; nwsteps = 0; break; } - if ( dns02[ii][jj] > 0 ) { nw = dns02[ii][jj]; nwsteps = 1.4142 * k; break; } - } - - - /* Southwest */ - for ( k = 1; k < infinit; k++ ) - { - ii = i + k; jj = j - k; - if (( ii == lines ) || ( jj == -1 )) { sw = 0; swsteps = 0; break; } - if ( dns02[ii][jj] > 0 ) { sw = dns02[ii][jj]; swsteps = 1.4142 * k; break; } - } - - - elev = 0; range = 0; - if ( nsteps > 0.5 ) { elev += north / nsteps; range += 1.00 / nsteps; } - if ( ssteps > 0.5 ) { elev += south / ssteps; range += 1.00 / ssteps; } - if ( esteps > 0.5 ) { elev += east / esteps; range += 1.00 / esteps; } - if ( wsteps > 0.5 ) { elev += west / wsteps; range += 1.00 / wsteps; } - if ( nwsteps > 0.5 ) { elev += nw / nwsteps; range += 1.00 / nwsteps; } - if ( nesteps > 0.5 ) { elev += ne / nesteps; range += 1.00 / nesteps; } - if ( swsteps > 0.5 ) { elev += sw / swsteps; range += 1.00 / swsteps; } - if ( sesteps > 0.5 ) { elev += se / sesteps; range += 1.00 / sesteps; } - - - temp = ( elev / range ) + 0.5 ; - if ( temp > 10000 ) temp = 10000; - if ( temp < 0 ) temp = 0; - - dns01[i][j] = temp ; - - } - } - - -/***************** CHECK FOR REMAINING ZEROES *********************/ - - -there: - -for ( i = 0; i < lines; i++ ) -for ( j = 0; j < samps; j++ ) -dns01 [i][j] = dns02 [i][j]; - - -if ( count > 3000 ) goto there2; /* Optional within-program early termination */ - - - for ( i = 0; i < lines; i++ ) - for ( j = 0; j < samps; j++ ) - if ( dns01 [i][j] == 0 ) goto again2; - - -if ( smooth == 1 ) goto there2; - - -/************************ SMOOTH THE RESULT ***********************/ - - - - dns03 = (unsigned short **) malloc( lines * sizeof( short * ) ); - if (dns03 == NULL) { printf("Not enough RAM.\a\n"); goto end; } - - for ( i = 0; i < lines; i ++ ) /* Get output RAM */ - { - dns03[ i ] = ( unsigned short * ) malloc( samps * sizeof(short) ); - if (dns03[ i ] == NULL) { printf("Not enough RAM.\a\n"); goto end; } - } - - - fseek( fp1, 0, SEEK_SET ); /* Reset pointer to Input Data */ - - /* Load original data into RAM Array 3 (Need data / Non-Data Info) */ - - /*vvvvvvvvvvvv Start Buffering of Input vvvvvvvvvvvv*/ - - buffy2 = (unsigned short *) malloc( samps*100 * sizeof(short) ); - if (buffy2 == NULL) { printf("Not enough RAM.\a\n"); goto end; } - - cycles = lines / 100; - for ( k = 0; k < cycles; k++ ) - { - fread( buffy2, sizeof(short), 100*samps, fp1 ); - for ( i = 0; i < 100; i++ ) - { - prevlines = ( k * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - dns03 [ prevlines ][ j ] = buffy2 [ prevpix + j ]; - } - } - - leftovers = lines - ( cycles * 100 ); - fread( buffy2, sizeof(short), samps * leftovers, fp1 ); - for ( i = 0; i < leftovers; i++ ) - { - prevlines = ( cycles * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - dns03 [ prevlines ][ j ] = buffy2 [ prevpix + j ]; - } - - /*^^^^^^^^^^^^^ End Buffering of Input ^^^^^^^^^^^^^*/ - - fprintf(stdout, "\n\nSmoothing Line Counter ( %5ld total ): ", lines ); - - - for ( i = 0; i < lines; i++ ) - { - for ( j = 0; j < samps; j++ ) - { - if ( dns03 [i][j] > 0 ) { dns02 [i][j] = dns01 [i][j]; continue; } - - - - for ( k = 1; k < infinit; k++ ) /* Find the smallest box size with data */ - { - for ( ii = i-k; ii <= i+k; ii++ ) - { - if ( ii < 0 ) continue; - if ( ii >= lines ) break; - for ( jj = j-k; jj <= j+k; jj++ ) - { - if ( jj < 0 ) continue; - if ( jj >= samps ) continue; - if ( dns03[ii][jj] > 0 ) goto k_nowset; - } - } - } - - k_nowset: - k = k / 4; if ( k < 1 ) k = 1; /* Errrrr. Make it fourth size */ - sum = 0; - for ( ii = i-k; ii <= i+k; ii++ ) - { - if ( ii < 0 ) iii = -1 - ii; else iii = ii; - if ( ii >= lines ) iii = lines - (ii - lines + 1); - for ( jj = j-k; jj <= j+k; jj++ ) - { - if ( jj < 0 ) jjj = -1 - jj; else jjj = jj; - if ( jj >= samps ) jjj = samps - (jj - samps + 1); - sum += dns01[iii][jjj]; - } - } - k = (k*2)+1; - dns02[i][j] = ( sum / k / k ) + 0.5; - - - - } - if ( !(i%100) ) - fprintf(stdout, "\b\b\b\b\b\b%5ld ", i ); /* Count lines on output screen */ - fflush(stdout); - } - fprintf(stdout, "\b\b\b\b\b\b%5ld ", lines ); - fflush(stdout); - - - - - - - -there2: - - - /*vvvvvvvvvvvv Start Buffering of Output vvvvvvvvvvvv*/ - - cycles = lines / 100; - for ( k = 0; k < cycles; k++ ) - { - for ( i = 0; i < 100; i++ ) - { - prevlines = ( k * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - buffy2 [ prevpix + j ] = dns02 [ prevlines ][ j ]; - } - fwrite( buffy2, sizeof(short), 100*samps, fp2 ); - } - - leftovers = lines - ( cycles * 100 ); - for ( i = 0; i < leftovers; i++ ) - { - prevlines = ( cycles * 100 ) + i; - prevpix = i * samps; - for ( j = 0; j < samps; j++ ) - buffy2 [ prevpix + j ] = dns02 [ prevlines ][ j ]; - } - fwrite( buffy2, sizeof(short), samps * leftovers, fp2 ); - - - /*^^^^^^^^^^^^^ End Buffering of Output ^^^^^^^^^^^^^*/ - - - fclose( fp1 ); fclose( fp2 ); - - time(&t2); printf ( "\n\nEnd "); printf( ctime(&t2) ); - printf( "\n Total time: "); - printf( "%4.2f minutes. \n", difftime(t2,t1)/60 ); - - - end: - printf( "Done.\a\n" ); - - return(0); - -} - -/*==================================================================================*/ - -int number_input() - -{ - char buffer[16]; - char *p; - - p=buffer; - while ( (*p = getc(stdin) ) != '\n') p++; - *p = '\0'; - return( atoi (buffer) ); -} - - - -/*==================================================================================*/ - -int yes_no_input() - -{ - char buffer[16]; - char *p; - - p=buffer; - while ( (*p = getc(stdin) ) != '\n') p++; - *p = '\0'; - if ( strlen(buffer) == 0 ) return ( 0 ); - return( 1 ); -} - - - - -/* EDITS: */ -/* */ -/* 03.04.02 Modified from Zero Fill Cycle. */ - Index: /issm/trunk/src/c/modules/HoleFillerx/HoleFillerx.cpp =================================================================== --- /issm/trunk/src/c/modules/HoleFillerx/HoleFillerx.cpp (revision 9312) +++ /issm/trunk/src/c/modules/HoleFillerx/HoleFillerx.cpp (revision 9313) @@ -1,6 +1,5 @@ /*!\file: HoleFillerx.cpp * \brief plug missing values of a double* array, using interpolation techniques. This - * is an adaptation of Dr. Crippen's routine (Jet Propulsion Lab, see HoleFillerCrippen.c - * for more details. + * is an adaptation of Dr. Crippen's routine (Jet Propulsion Lab) */