---

median.c

Go to the documentation of this file.

#include        "sadie.h"

/*-Copyright Information------------------------------------------------------*/
/* Copyright (c) 1988 by the University of Arizona Digital Image Analysis Lab */
/*----------------------------------------------------------------------------*/
/*-General Information--------------------------------------------------------*/
/*                                                                            */
/*   This function performs median filtering of an image.                     */
/*                                                                            */
/*----------------------------------------------------------------------------*/
/*-Interface Information------------------------------------------------------*/
void MEDIAN (
IMAGE *in,      /*  I   Pointer to the input image.                           */
short jsize,    /*  I   Number of lines in the median window.                 */
short ksize,    /*  I   Number of pixels/line in the median window.           */
IMAGE **out     /*  O   Address of a pointer to the output image.             */
/*----------------------------------------------------------------------------*/
) { register short i, j, k, l, m, count, index, jstart, jend, kstart, kend, jover=(jsize-1)/2, kover=(ksize-1)/2;
    char   msg[SLEN];
    double pinc, psum;
    PIXEL  med, *list=0;

    if (TIMES) TIMING(T_MEDIAN);
    if (NAMES) {
        MESSAGE('I',"");
        MESSAGE('I',"MEDIAN");
        MESSAGE('I',"");
        sprintf(msg," Input image:           %s",in->text);
        MESSAGE('I',msg);
        sprintf(msg," Window size: lines:    %d",jsize);
        MESSAGE('I',msg);
        sprintf(msg,"              pixels:   %d",ksize);
        MESSAGE('I',msg);
        MESSAGE('I'," ...............");
    }

    /* check input */
    if (!CHECKIMG(in)) {
        MESSAGE('E'," Can't identify image.");
        goto the_end;
    } else if (jsize <= 0  ||  jsize%2 == 0) {
        MESSAGE('E'," Number of lines in the median window must be an odd number greater than zero.");
        goto the_end;
    } else if (jsize > in->nlin) {
        MESSAGE('E'," Image line size must be equal to or greater than window size.");
        goto the_end;
    } else if (ksize <= 0  ||  ksize%2 == 0) {
        MESSAGE('E'," Number of pixels/line in the median window must be an odd number greater than zero.");
        goto the_end;
    } else if (ksize > in->npix) {
        MESSAGE('E'," Image pixel size must be equal to or greater than window size.");
        goto the_end;
    }

    /* create image of appropriate size */
    if (!CHECKIMG(*out)) GETMEM(in->nbnd,in->nlin,in->npix,out);
    if (!*out) goto the_end;

    /* allocate space for 1-D array */
    if (!(list=(PIXEL *)malloc(jsize*ksize*sizeof(PIXEL)))) {
        MESSAGE('E'," Memory request failed.");
        goto the_end;
    }

    /* initialize progress indicator */
    if (LINES  &&  !PROGRESS(psum=0.0)) goto the_end;
    pinc = 1.0/(double)in->nbnd/(double)in->nlin;

    /* put pixels inside window into array */
    for (i=0; i<in->nbnd; i++) {
        for (j=0; j<in->nlin; j++) {
            jstart = max(0,j-jover);
            jend   = min(j+jover,in->nlin-1);
            for (k=0; k<in->npix; k++) {
                kstart = max(0,k-kover);
                kend   = min(k+kover,in->npix-1);
                for (count=0,l=jstart; l<=jend; l++) {
                    for (m=kstart; m<=kend; m++) {
                        list[count++] = in->data[i][l][m];
                    }
                }

                /* compute median */
                for (l=0; l<(count+1)/2; l++) {
                    med = list[index=l];
                    for (m=l+1; m<count; m++) {
                        if (med > list[m]) med = list[index=m];
                    }
                    list[index] = list[l];
                }
                (*out)->data[i][j][k] = med;
            }
            if (LINES  &&  !PROGRESS(psum+=pinc)) goto the_end;
        }
    }

    the_end:
    if (list)  free(list);
    if (TIMES) TIMING(T_EXIT);
}

---