Main Page   Data Structures   File List   Data Fields   Globals  

ROIsconvl.c

Go to the documentation of this file.
00001 #include        "sadie.h"
00002 #include        "math.h"
00003 
00004 
00005 typedef struct {
00006     unsigned short hstart;
00007     unsigned short hend;
00008     unsigned short *vstart;
00009     unsigned short *vend;
00010 }                    LOCAL_INDEX;
00011 
00012 
00013 /*-Copyright Information------------------------------------------------------*/
00014 /* Copyright (c) 1988 by the University of Arizona Digital Image Analysis Lab */
00015 /*----------------------------------------------------------------------------*/
00016 /*-General Information--------------------------------------------------------*/
00017 /*                                                                            */
00018 /*   This function computes a two-dimensional spatial convolution between     */
00019 /*   a user-supplied spatial kernel and a single-band image.                  */
00020 /*                                                                            */
00021 /*----------------------------------------------------------------------------*/
00022 /*-Interface Information------------------------------------------------------*/
00023 void ROI_SCONVL (
00024 IMAGE *in,      /*  I   Pointer to the input image.                           */
00025 LOCAL_INDEX *localindex,                 
00026 PIXEL *psf,     /*  I   Pointer to the convolution mask matrix[jsize][ksize]. */
00027 short jsize,    /*  I   Number of lines in the convolution mask.              */
00028 short ksize,    /*  I   Number of pixels/line in the convolution mask.        */
00029 IMAGE **out     /*  O   Address of a pointer to the output image.             */
00030 /*----------------------------------------------------------------------------*/
00031 ) { register int j, k, l, m, n, jmarg=jsize/2, kmarg=ksize/2, box;
00032     char   msg[SLEN];
00033     double sum, r1, rn, pinc, psum;
00034 
00035     int progmeter;
00036     double progress;
00037     int progcount, progtotal, progincr;
00038     int cancel=0;
00039 
00040     progmeter = ProgMeter_Create("Convolution");
00041     progress = 0.0;
00042 
00043 
00044     if (TIMES) TIMING(T_SCONVL);
00045     if (NAMES) {
00046         MESSAGE('I',"");
00047         MESSAGE('I',"SCONVL");
00048         MESSAGE('I',"");
00049         sprintf(msg," Input image:           %s",in->text);
00050         MESSAGE('I',msg);
00051         sprintf(msg," Kernel size: lines:    %d",jsize);
00052         MESSAGE('I',msg);
00053         sprintf(msg,"              pixels:   %d",ksize);
00054         MESSAGE('I',msg);
00055         if (jsize <= SLEN/16  &&  ksize <= SLEN/16) {
00056             MESSAGE('I'," Spatial kernel:");
00057             for (n=j=0; j<jsize; j++) {
00058                 for (m=k=0,msg[m++]=' '; k<ksize; k++,m=strlen(msg)) {
00059                     sprintf(msg+m,"%12.4e",(double)psf[n++]);
00060                 }
00061                 MESSAGE('I',msg);
00062             }
00063         } else {
00064             MESSAGE('I'," Kernel too large to list.");
00065         }
00066         MESSAGE('I'," ...............");
00067     }
00068 
00069     /* check input */
00070     if (!CHECKIMG(in)) {
00071         MESSAGE('E'," Can't identify image.");
00072         goto the_end;
00073     } else if (in->nbnd > 1) {
00074         MESSAGE('E'," Image must be single-band.");
00075         goto the_end;
00076     } else if (jsize < 0  ||  jsize%2 == 0) {
00077         MESSAGE('E'," Number of lines in the convolution mask must be an odd number greater than zero.");
00078         goto the_end;
00079     } else if (jsize > in->nlin) {
00080         MESSAGE('E'," Image size must be equal to or greater than convolution mask size in the line dimension.");
00081         goto the_end;
00082     } else if (ksize < 0  ||  ksize%2 == 0) {
00083         MESSAGE('E'," Number of pixels/line in the convolution mask must be an odd number greater than zero.");
00084         goto the_end;
00085     } else if (ksize > in->npix) {
00086         MESSAGE('E'," Image size must be equal to or greater than convolution mask size in the pixel dimension.");
00087         goto the_end;
00088     }
00089 
00090     /* create image of appropriate size */
00091     if (!CHECKIMG(*out)) GETMEM(1,in->nlin,in->npix,out);
00092     if (!*out) goto the_end;
00093 
00094     /* Initialize output image */
00095     for(j = 0; j < in->nlin; j++) {
00096         for(k = 0; k < in->npix; k++) {
00097             (*out)->data[0][j][k] = (PIXEL)0;
00098         }
00099     }
00100 
00101     /* initialize progress indicator */
00102     if (LINES  &&  !PROGRESS(psum=0.0)) goto the_end;
00103     pinc = 1.0/(double)(in->nlin-2*jmarg);
00104     progtotal = in->nlin-2*jmarg;
00105     progcount = 0;
00106     progincr = rint((double)progtotal/20.0);
00107 
00108     /* see if it's a box-filter PSF */
00109     for (box=TRUE,k=1; box && k<jsize*ksize; k++) {
00110         box = psf[k] == psf[0];
00111     }
00112     
00113     if (box) { /* Box filter implementation is yet to be tested */
00114         for (k = localindex->hstart + kmarg; k <= localindex->hend - kmarg; k++) {
00115 
00116             /* initialize box */
00117             for (r1 = sum  = 0.0, l = k -kmarg; l <= k + kmarg; l++) {
00118                 r1 += (double)in->data[0][localindex->vstart[k] - jmarg][l];
00119                 for (j = localindex->vstart[k] - jmarg + 1; j <= localindex->vstart[k] + jmarg; j++) {
00120                     sum += (double)in->data[0][j][l];
00121                 }
00122             }
00123             (*out)->data[0][localindex->vstart[k]][k] = (PIXEL)((sum += r1)*(double)psf[0]);
00124 
00125             /* box-filter algorithm */
00126             for (j = localindex->vstart[k] + 1; j <= localindex->vend[k]; j++) {
00127                 for (sum -= r1, r1 = rn = 0.0, l = k - kmarg; l <= k + kmarg; l++) {
00128                     r1 += (double)in->data[0][j - jmarg][l];
00129                     rn += (double)in->data[0][j + jmarg][l];
00130                 }
00131                 (*out)->data[0][j][k] = (PIXEL)((sum += rn)*(double)psf[0]);
00132             }
00133             if (LINES  &&  !PROGRESS(psum+=pinc)) goto the_end;
00134             if (progmeter) {
00135                     progcount++;
00136                     if ((progcount%progincr) == 0) {
00137                         progress = ((double)progcount/(double)progtotal)*100.0;
00138                         cancel = ProgMeter_Update(progmeter,progress);
00139                         if (cancel != 0) goto the_end;
00140                     }
00141             }
00142         }
00143 
00144     } else /* compute normal convolution */ {
00145         for (k = localindex->hstart; k <= localindex->hend; k++) {
00146             for (j = localindex->vstart[k]; j <= localindex->vend[k]; j++) {
00147                 for (sum = 0.0,n = jsize*ksize, l = 0; l < jsize; l++) {
00148                     for (m = 0; m < ksize; m++) {
00149                         sum += (double)in->data[0][j-jmarg+l][k-kmarg+m]*(double)psf[--n];
00150                     }
00151                 }
00152                 (*out)->data[0][j][k] = (PIXEL)sum;
00153             }
00154             if (LINES  &&  !PROGRESS(psum+=pinc)) goto the_end;
00155             if (progmeter) {
00156                     progcount++;
00157                     if ((progcount%progincr) == 0) {
00158                         progress = ((double)progcount/(double)progtotal)*100.0;
00159                         cancel = ProgMeter_Update(progmeter,progress);
00160                         if (cancel != 0) goto the_end;
00161                     }
00162             }
00163         }
00164     }
00165 
00166     
00167     the_end:
00168     if (TIMES) TIMING(T_EXIT);
00169 
00170     if (progmeter) {
00171         ProgMeter_Destroy(progmeter);
00172     }
00173 
00174     if (cancel != 0) {
00175         if (*out) RELMEM(*out);
00176         *out = NULL;
00177     }
00178 }

Generated on Wed Apr 9 08:56:10 2003 for TREES by doxygen1.2.14 written by Dimitri van Heesch, © 1997-2002