find_register.c Source File

00001 #include "sadie.h"
00002 #include "proto.h"
00003 #include <math.h>
00004 #include <sys/time.h>
00005 
00006 
00007 /*----------------------------------------------------------------------------*/
00008 /*-General Information--------------------------------------------------------*/
00009 /*                                                                            */
00010 /* This function performs spatial correlations between two images in order    */
00011 /* to determine rough registration offsets.  The input images are averaged    */
00012 /* and subsampled by a "shrink" factor in order to decrease the sensitivity   */
00013 /* to errors in the estimated x and y offsets.                                */
00014 /*                                                                            */
00015 /* This procedure allows for an error of up to                                */
00016 /*       (search_size - target_size)*shrink/2 pixels                          */
00017 /* in the estimated x and y offsets.                                          */
00018 /*                                                                            */
00019 /*----------------------------------------------------------------------------*/
00020 /*-Interface Information------------------------------------------------------*/
00021 void FIND_REGISTER_ROUGH (
00022 IMAGE  *img1,       /*  I   Pointer to the first input image.                    */
00023 IMAGE  *img2,       /*  I   Pointer to the second input image.                   */
00024 int    est_x_off,   /*  I   Estimated x-offset of second image                   */
00025 int    est_y_off,   /*  I   Estimated y-offset of second image                   */
00026 int    search_size, /*  I   Search window size in first image                    */
00027 int    target_size, /*  I   Target window size in second image                   */
00028 int    shrink,      /*  I   Factor by which to shrink images                     */
00029 int    *x_off,      /*  O   Actual x-offset of second image                      */
00030 int    *y_off       /*  O   Actual y-offset of second image                      */
00031 /*-------------------------------------------------------------------------------*/
00032 ) { register int  i, j, k, m, n;
00033     char msg[SLEN];
00034     int x_slop, y_slop;
00035     int roi1_x1, roi1_x2, roi1_y1, roi1_y2;
00036     int roi2_x1, roi2_x2, roi2_y1, roi2_y2, roi2_x_size, roi2_y_size;
00037     int targ_x_off, targ_y_off;
00038     int npixval, max_i, max_j;
00039     double df, mean_t, old_mean, var_t, sd_t, sd_s;
00040     double sum_s,         ss_s;
00041     double sum_lap,       ss_lap;
00042     double sum_oldline,   ss_oldline;
00043     double sum_lastline,  ss_lastline;
00044     double sum_firstline, ss_firstline;
00045     double oldval, newval;
00046     double r, sum_prod;
00047     double pixval, maxval;
00048     IMAGE *small_inimg[2];
00049 
00050     x_slop = y_slop = (search_size - target_size);
00051     for (i=0; i<2; i++) {
00052         small_inimg[i] = NULL;
00053     }
00054 
00055     /* check input */
00056     if (!CHECKIMG(img1)) {
00057         MESSAGE('E'," Can't identify first input image.");
00058         goto the_end;
00059         }
00060     if (!CHECKIMG(img2)) {
00061         MESSAGE('E'," Can't identify second input image.");
00062         goto the_end;
00063         }
00064     
00065     /* Shrink input images */
00066     RESAMPL(img1,shrink,shrink,shrink,shrink,&small_inimg[0]);
00067     RESAMPL(img2,shrink,shrink,shrink,shrink,&small_inimg[1]);
00068     if (!CHECKIMG(small_inimg[0])) {
00069         MESSAGE('E'," Can't identify first small input image.");
00070         goto the_end;
00071         }
00072     if (!CHECKIMG(small_inimg[1])) {
00073         MESSAGE('E'," Can't identify second small input image.");
00074         goto the_end;
00075         }
00076     
00077     est_x_off = rint((double)est_x_off/(double)shrink);
00078     est_y_off = rint((double)est_y_off/(double)shrink);
00079 
00080     if (est_x_off > small_inimg[0]->npix) {
00081             MESSAGE('E'," Images cannot be registered to one another because they do not overlap!");
00082             goto the_end;
00083     }
00084     if (est_y_off > small_inimg[0]->nlin) {
00085             MESSAGE('E'," Images cannot be registered to one another because they do not overlap!");
00086             goto the_end;
00087     }
00088 
00089     /* Determine bounds on estimated overlapping regions */
00090     if (est_x_off >= 0) {
00091             roi1_x1 = min(est_x_off,small_inimg[0]->npix-1);
00092             roi1_x2 = min(small_inimg[0]->npix-1,(small_inimg[1]->npix-1)+roi1_x1);
00093 
00094             roi2_x1 = 0;
00095             roi2_x2 = min(small_inimg[1]->npix-1,(small_inimg[0]->npix-1)-est_x_off);
00096     } else {
00097             roi1_x1 = 0;
00098             roi1_x2 = min(small_inimg[0]->npix-1,(small_inimg[1]->npix-1) + est_x_off);
00099 
00100             roi2_x1 = min(small_inimg[1]->npix-1, 0-est_x_off);
00101             roi2_x2 = min(small_inimg[1]->npix-1, (small_inimg[0]->npix-1) - est_x_off);
00102     }
00103 
00104     if (est_y_off >= 0) {
00105             roi1_y1 = min(est_y_off,small_inimg[0]->nlin-1);
00106             roi1_y2 = min(small_inimg[0]->nlin-1,(small_inimg[1]->nlin-1)+roi1_y1);
00107 
00108             roi2_y1 = 0;
00109             roi2_y2 = min(small_inimg[1]->nlin-1, (small_inimg[0]->nlin-1)-est_y_off);
00110     } else {
00111             roi1_y1 = 0;
00112             roi1_y2 = min(small_inimg[0]->nlin-1,(small_inimg[1]->nlin-1) + est_y_off);
00113 
00114             roi2_y1 = min(small_inimg[1]->nlin-1, 0-est_y_off);
00115             roi2_y2 = min(small_inimg[1]->nlin-1, (small_inimg[0]->nlin-1) - est_y_off);
00116     }   
00117 
00118     /* Make sure there is enough overlapping data to work with */
00119     if (((roi1_x2 - roi1_x1) < search_size)||((roi1_y2 - roi1_y1) < search_size)) {
00120             MESSAGE('E'," Images cannot be registered to one another because the overlap area is too small!");
00121             goto the_end;
00122     }
00123     if (((roi1_x2 - roi1_x1) != (roi2_x2 - roi2_x1))||((roi1_y2 - roi1_y1) != (roi2_y2 - roi2_y1))) {
00124             MESSAGE('E'," Overlapping areas are not the same size!");
00125             goto the_end;
00126     }
00127 
00128     /* Compute Maximum Cross-Correlation Coefficient */
00129 
00130     /* Re-size the ROI in the target to give multiple potential matches */
00131     targ_x_off = x_slop / 2;
00132     targ_y_off = y_slop / 2;
00133     roi2_x1 += targ_x_off;
00134     roi2_x2 -= (x_slop - targ_x_off);
00135     roi2_x_size = 1 + roi2_x2 - roi2_x1;
00136     roi2_y1 += targ_y_off;
00137     roi2_y2 -= (y_slop - targ_y_off);
00138     roi2_y_size = 1 + roi2_y2 - roi2_y1;
00139 
00140     printf("roi1 = (%d,%d) to (%d,%d)\n",roi1_x1,roi1_y1,roi1_x2,roi1_y2);
00141     printf("roi2 = (%d,%d) to (%d,%d)\n",roi2_x1,roi2_y1,roi2_x2,roi2_y2);
00142 
00143     /* Compute target size, mean and standard deviation */
00144     var_t = 0;
00145     n = roi2_x1;
00146     mean_t = small_inimg[1]->data[0][roi2_y1][n++];
00147     npixval = 1;
00148     for (m=roi2_y1; m<=roi2_y2; m++) {
00149         while (n <= roi2_x2) {
00150             pixval = small_inimg[1]->data[0][m][n++];
00151             old_mean = mean_t;
00152             mean_t = old_mean + (pixval - old_mean)/(++npixval);
00153             var_t += (pixval - old_mean)*(pixval - mean_t);
00154         }
00155         n = roi2_x1;
00156     }
00157     df = npixval - 1;
00158     sd_t = sqrt(var_t / df);
00159     /* initialize sums for first search position */
00160     sum_lap = 0;
00161     ss_lap = 0;
00162     for (m=roi1_y1; m<(roi1_y1+roi2_y_size); m++) {
00163         for (n=roi1_x1; n<(roi1_x1+roi2_x_size); n++) {
00164             pixval = small_inimg[0]->data[0][m][n];
00165             sum_lap += pixval;
00166             ss_lap  += pixval * pixval;
00167         }
00168     }
00169     sum_oldline = 0;
00170     ss_oldline  = 0;
00171     maxval = -1;
00172     max_i = max_j = 0;
00173     for (j=0; j<y_slop; j++) {
00174         sum_lastline = 0;   /* sums for line added to end of search area */
00175         ss_lastline  = 0;
00176         sum_firstline = 0;  /* sums for line next to drop from search area */
00177         ss_firstline  = 0;
00178         for (k=roi1_x1; k<(roi1_x1+roi2_x_size); k++) {
00179             pixval = small_inimg[0]->data[0][roi1_y1+j+roi2_y_size-1][k];
00180             sum_lastline += pixval;
00181             ss_lastline  += pixval * pixval;
00182             pixval = small_inimg[0]->data[0][roi1_y1+j][k];
00183             sum_firstline += pixval;
00184             ss_firstline  += pixval * pixval;
00185         }
00186         sum_s = sum_lap; /* initialize sum for this y position */ 
00187         ss_s  = ss_lap;
00188         sum_lap += sum_lastline - sum_oldline; /* update for next y positon */
00189         ss_lap  += ss_lastline  - ss_oldline;
00190         sum_oldline = sum_firstline;
00191         ss_oldline  = ss_firstline;
00192         /* searching first x position treated as a special case */
00193         sum_prod = 0;
00194         for (m=0; m<roi2_y_size; m++) {
00195             for (n=0; n<roi2_x_size; n++) {
00196                 sum_prod += small_inimg[0]->data[0][roi1_y1+j+m][roi1_x1+n]
00197                           * small_inimg[1]->data[0][roi2_y1+m][roi2_x1+n];
00198             }
00199         }
00200         sum_prod -= sum_s*mean_t;
00201         r = sum_prod / (df * sd_t * sqrt((npixval * ss_s - sum_s*sum_s)/(npixval*df)));
00202         if (r > maxval) {
00203             maxval = r;
00204             max_i = 0;
00205             max_j = j;
00206         }
00207         /* searching subsequent x positions */
00208         for (i=1; i<x_slop; i++) {
00209             sum_prod = 0;
00210             for (m=0; m<roi2_y_size; m++) {
00211                 newval = small_inimg[0]->data[0][roi1_y1+j+m][roi1_x1+i+roi2_x_size-1];
00212                 oldval = small_inimg[0]->data[0][roi1_y1+j+m][roi1_x1+i-1];
00213                 sum_s += (newval - oldval);
00214                 ss_s  += (newval*newval - oldval*oldval);
00215                 for (n=0; n<roi2_x_size; n++) {
00216                     sum_prod += small_inimg[0]->data[0][roi1_y1+j+m][roi1_x1+i+n]
00217                               * small_inimg[1]->data[0][roi2_y1+m][roi2_x1+n];
00218                 }
00219             }
00220             sd_s = sqrt((npixval * ss_s - sum_s*sum_s)/(npixval*df));
00221             sum_prod -= sum_s*mean_t;
00222             r = sum_prod / (df * sd_t * sd_s);
00223             // printf("r=%8.5f, sd_t=%8.5f, sd_s=%8.5f\n", r, sd_t, sd_s);
00224             if (r > maxval) {
00225                 maxval = r;
00226                 max_i = i;
00227                 max_j = j;
00228             }
00229         }
00230     }
00231 
00232     *x_off = shrink * (est_x_off + max_i - targ_x_off);
00233     *y_off = shrink * (est_y_off + max_j - targ_y_off);
00234 
00235     the_end:
00236     
00237     if (CHECKIMG(small_inimg[0])) {
00238         RELMEM(small_inimg[0]);
00239     }
00240 
00241     if (CHECKIMG(small_inimg[1])) {
00242         RELMEM(small_inimg[1]);
00243     }
00244 
00245 }
00246 
00247 
00248 
00249 
00250 /*----------------------------------------------------------------------------*/
00251 /*-General Information--------------------------------------------------------*/
00252 /*                                                                            */
00253 /* This function performs spatial correlations between two images in order    */
00254 /* to determine the best registration offsets.  To improve robustness,        */
00255 /* the target and search areas are split into four separate regions.          */
00256 /* See W. S. Conner's thesis for details.                                     */
00257 /*                                                                            */
00258 /*----------------------------------------------------------------------------*/
00259 /*-Interface Information------------------------------------------------------*/
00260 void FIND_REGISTER (
00261 IMAGE  *img1,     /*  I   Pointer to the first input image.                    */
00262 IMAGE  *img2,     /*  I   Pointer to the second input image.                   */
00263 int    est_x_off, /*  I   Estimated x-offset of second image                   */
00264 int    est_y_off, /*  I   Estimated y-offset of second image                   */
00265 int    *x_off,    /*  I   Actual x-offset of second image                      */
00266 int    *y_off,    /*  I   Actual y-offset of second image                      */
00267 PIXEL  *bias_adj, /*  I   Bias Adjustment of second image                      */
00268 PIXEL  *gain_adj  /*  I   Gain Adjustment of second image                      */
00269 /*----------------------------------------------------------------------------*/
00270 ) { register int i, j, k, m, n;
00271     int max_i, max_j;
00272     char msg[SLEN];
00273     int x_index, y_index;
00274     int x_offset[4], y_offset[4];
00275     int roi1_x1, roi1_x2, roi1_y1, roi1_y2;
00276     int roi2_x1, roi2_x2, roi2_y1, roi2_y2;
00277     int roi1_center_x, roi1_center_y, roi2_center_x, roi2_center_y;
00278     PIXEL mean_s, mean_t;
00279     PIXEL denom_s, denom_t, num;
00280     PIXEL temp_s, temp_t, value, maxval;
00281 
00282     /*
00283      * This procedure allows for an error in estimated offsets of up to
00284      *     (search_size - target_size) / 2
00285      * pixels in both the x and y directions.
00286      */
00287     int search_size = 49;
00288     int target_size = 29;
00289     
00290     int rough_search_size = 69;
00291     int rough_target_size = 29;
00292     int rough_shrink = 5;
00293 
00294     IMAGE *inimg[2];
00295     int jbgn[2], kbgn[2];
00296     int nlin, npix;
00297     int xcount, ycount, counter;
00298     struct  timeval start, find_first_reg, findreg, findgain, end;
00299     int est2_x_off = 0;
00300     int est2_y_off = 0;
00301         
00302     /* start timing */
00303     gettimeofday(&start, NULL);
00304 
00305     /* check input */
00306     if (!CHECKIMG(img1)) {
00307         MESSAGE('E'," Can't identify first input image.");
00308         goto the_end;
00309         }
00310     if (!CHECKIMG(img2)) {
00311         MESSAGE('E'," Can't identify second input image.");
00312         goto the_end;
00313         }
00314     
00315     inimg[0] = img1;
00316     inimg[1] = img2;
00317         
00318 
00319     if (NAMES) {
00320         MESSAGE('I',"");
00321         MESSAGE('I',"REGISTER");
00322         MESSAGE('I',"");
00323         sprintf(msg," First input image:                   %s",img1->text);
00324         MESSAGE('I',msg);
00325         sprintf(msg," Second input image:                  %s",img2->text);
00326         MESSAGE('I',msg);
00327         MESSAGE('I',"");
00328         MESSAGE('I',"");
00329         sprintf(msg," Estimated x-offset:  %d",est_x_off);
00330                 MESSAGE('I',msg);
00331         sprintf(msg," Estimated y-offset:  %d",est_y_off);
00332                 MESSAGE('I',msg);
00333     }
00334 
00335     
00336     if (est_x_off > img1->npix) {
00337             MESSAGE('E'," Images cannot be registered to one another because they do not overlap!");
00338             goto the_end;
00339     }
00340     if (est_y_off > img1->nlin) {
00341             MESSAGE('E'," Images cannot be registered to one another because they do not overlap!");
00342             goto the_end;
00343     }
00344 
00345 
00346     
00347 
00348 
00349 
00350     /* note time at beginning of first correlation */
00351     gettimeofday(&find_first_reg, NULL);
00352 
00353 
00354     /*
00355      * Refine the estimated x and y offsets by bringining them to within
00356      * rough_shrink pixels through an initial rough registration.
00357      */
00358     FIND_REGISTER_ROUGH (img1,img2,est_x_off,est_y_off,rough_search_size,rough_target_size,rough_shrink,&est2_x_off,&est2_y_off);
00359     //IMG2DISK(img1,32,"img1.SAD");
00360     //IMG2DISK(img2,32,"img2.SAD");
00361 
00362 printf("Rough Registration results: (%d,%d) ==> (%d,%d)\n",est_x_off, est_y_off, est2_x_off, est2_y_off);
00363 
00364     if (est2_x_off > img1->npix) {
00365             MESSAGE('E'," Images cannot be registered to one another because they do not overlap!");
00366             goto the_end;
00367     }
00368     if (est2_y_off > img1->nlin) {
00369             MESSAGE('E'," Images cannot be registered to one another because they do not overlap!");
00370             goto the_end;
00371     }
00372 
00373     /* Determine bounds on estimated overlapping regions */
00374     if (est2_x_off >= 0) {
00375             roi1_x1 = min(est2_x_off,img1->npix-1);
00376             roi1_x2 = min(img1->npix-1,(img2->npix-1)+roi1_x1);
00377 
00378             roi2_x1 = 0;
00379             roi2_x2 = min(img2->npix-1,(img1->npix-1)-est2_x_off);
00380     } else {
00381             roi1_x1 = 0;
00382             roi1_x2 = min(img1->npix-1,(img2->npix-1) + est2_x_off);
00383 
00384             roi2_x1 = min(img2->npix-1, 0-est2_x_off);
00385             roi2_x2 = min(img2->npix-1, (img1->npix-1) - est2_x_off);
00386     }
00387 
00388     if (est2_y_off >= 0) {
00389             roi1_y1 = min(est2_y_off,img1->nlin-1);
00390             roi1_y2 = min(img1->nlin-1,(img2->nlin-1)+roi1_y1);
00391 
00392             roi2_y1 = 0;
00393             roi2_y2 = min(img2->nlin-1, (img1->nlin-1)-est2_y_off);
00394     } else {
00395             roi1_y1 = 0;
00396             roi1_y2 = min(img1->nlin-1,(img2->nlin-1) + est2_y_off);
00397 
00398             roi2_y1 = min(img2->nlin-1, 0-est2_y_off);
00399             roi2_y2 = min(img2->nlin-1, (img1->nlin-1) - est2_y_off);
00400     }   
00401 
00402     /* Make sure there is enough overlapping data to work with */
00403     if (((roi1_x2 - roi1_x1) < (2*search_size))||((roi1_y2 - roi1_y1) < (2*search_size))) {
00404             MESSAGE('E'," Images cannot be registered to one another because the overlap area is too small!");
00405             goto the_end;
00406     }
00407     if (((roi1_x2 - roi1_x1) != (roi2_x2 - roi2_x1))||((roi1_y2 - roi1_y1) != (roi2_y2 - roi2_y1))) {
00408             MESSAGE('E'," Overlapping areas are not the same size!");
00409             goto the_end;
00410     }
00411 
00412 
00413     /* note time at beginning of main correlation */
00414     gettimeofday(&findreg, NULL);
00415 
00416     mean_t = (PIXEL)0.0;
00417     mean_s = (PIXEL)0.0;
00418 
00419     /* Split target and search areas into four separate regions */
00420     /*   See W.S. Conner's thesis for details.                  */
00421     for (xcount = counter = 0; xcount <2; xcount++) {
00422         for (ycount = 0; ycount < 2; ycount++,counter++) {
00423         
00424             roi1_center_x = roi1_x1 + ((roi1_x2 - roi1_x1)/4) + (xcount*((roi1_x2 - roi1_x1)/2));
00425             roi1_center_y = roi1_y1 + ((roi1_y2 - roi1_y1)/4) + (ycount*((roi1_y2 - roi1_y1)/2));
00426             roi2_center_x = roi2_x1 + ((roi2_x2 - roi2_x1)/4) + (xcount*((roi2_x2 - roi2_x1)/2));
00427             roi2_center_y = roi2_y1 + ((roi2_y2 - roi2_y1)/4) + (ycount*((roi2_y2 - roi2_y1)/2));
00428             
00429             
00430             /* Compute mean target value */
00431             for (m=0; m<target_size; m++) {
00432                     for (n=0; n<target_size; n++) {
00433                                 mean_t += img2->data[0][(roi2_center_y-((target_size-1)/2))+m][(roi2_center_x-((target_size-1)/2))+n];
00434                 }
00435             }
00436 
00437             /* Compute mean search value */
00438             for (m=0; m<search_size; m++) {
00439                     for (n=0; n<search_size; n++) {
00440                                 mean_s += img1->data[0][(roi1_center_y-((search_size-1)/2))+m][(roi1_center_x-((search_size-1)/2))+n];
00441                     }
00442             }
00443 
00444         }
00445     }
00446     mean_t /= (PIXEL)(4 * target_size * target_size);
00447     mean_s /= (PIXEL)(4 * search_size * search_size);
00448     
00449     denom_t = (PIXEL)0.0;
00450 
00451     /* Split target and search areas into four separate regions */
00452     /*   See W.S. Conner's thesis for details.                  */
00453     for (xcount = counter = 0; xcount <2; xcount++) {
00454         for (ycount = 0; ycount < 2; ycount++,counter++) {
00455         
00456             roi1_center_x = roi1_x1 + ((roi1_x2 - roi1_x1)/4) + (xcount*((roi1_x2 - roi1_x1)/2));
00457             roi1_center_y = roi1_y1 + ((roi1_y2 - roi1_y1)/4) + (ycount*((roi1_y2 - roi1_y1)/2));
00458             roi2_center_x = roi2_x1 + ((roi2_x2 - roi2_x1)/4) + (xcount*((roi2_x2 - roi2_x1)/2));
00459             roi2_center_y = roi2_y1 + ((roi2_y2 - roi2_y1)/4) + (ycount*((roi2_y2 - roi2_y1)/2));
00460             
00461             
00462             /* Compute the target denominator value */
00463             for (m=0; m<target_size; m++) {
00464                     for (n=0; n<target_size; n++) {
00465                                 denom_t += (img2->data[0][(roi2_center_y-((target_size-1)/2))+m][(roi2_center_x-((target_size-1)/2))+n] - mean_t)*(img2->data[0][(roi2_center_y-((target_size-1)/2))+m][(roi2_center_x-((target_size-1)/2))+n] - mean_t);
00466                     }
00467             }
00468 
00469         }
00470     }
00471     denom_t = (PIXEL) sqrt((double)denom_t);
00472 
00473 
00474     //printf("mean_t = %f\n",mean_t);
00475     //printf("mean_s = %f\n",mean_s);
00476     //printf("denom_t = %f\n",denom_t);
00477 
00478     maxval = (PIXEL)0.0;
00479     max_i = max_j = 0;
00480     for (j=0; j<(search_size - target_size)+1; j++) {
00481         for (i=0; i<(search_size - target_size)+1; i++) {
00482                 denom_s = (PIXEL)0.0;
00483             num = (PIXEL)0.0;
00484                 
00485             /* Split target and search areas into four separate regions */
00486             /*   See W.S. Conner's thesis for details.                  */
00487             for (xcount = counter = 0; xcount <2; xcount++) {
00488                 for (ycount = 0; ycount < 2; ycount++,counter++) {
00489                     roi1_center_x = roi1_x1 + ((roi1_x2 - roi1_x1)/4) + (xcount*((roi1_x2 - roi1_x1)/2));
00490                     roi1_center_y = roi1_y1 + ((roi1_y2 - roi1_y1)/4) + (ycount*((roi1_y2 - roi1_y1)/2));
00491                     roi2_center_x = roi2_x1 + ((roi2_x2 - roi2_x1)/4) + (xcount*((roi2_x2 - roi2_x1)/2));
00492                     roi2_center_y = roi2_y1 + ((roi2_y2 - roi2_y1)/4) + (ycount*((roi2_y2 - roi2_y1)/2));
00493 
00494                     for (m=0; m<target_size; m++) {
00495                                 for (n=0; n<target_size; n++) {
00496 
00497                                     temp_s = img1->data[0][(roi1_center_y-((search_size-1)/2))+m+j][(roi1_center_x-((search_size-1)/2))+n+i] - mean_s;
00498                                     temp_t = img2->data[0][(roi2_center_y-((target_size-1)/2))+m][(roi2_center_x-((target_size-1)/2))+n] - mean_t;
00499 
00500                                     denom_s += temp_s*temp_s;
00501                                     num += temp_t*temp_s;
00502                         }
00503                     }
00504                 }
00505                 }
00506             
00507             denom_s = (PIXEL) sqrt((double)denom_s);
00508 
00509             value = num / (denom_t * denom_s);
00510 
00511 
00512             if (value > maxval) {
00513                 maxval = value;
00514                 max_i = i;
00515                 max_j = j;
00516             }
00517 
00518         }
00519     }
00520 
00521     *x_off = est2_x_off + (max_i - ((search_size-target_size)/2));
00522     *y_off = est2_y_off + (max_j - ((search_size-target_size)/2));
00523 
00524     //printf("offset = (%d,%d)\n", *x_off, *y_off);
00525     //printf("max = (%d,%d)\n", max_i, max_j);
00526     fprintf(stdout,"deltamax = (%d,%d)\n", (max_i - ((search_size-target_size)/2)), (max_j - ((search_size-target_size)/2)));
00527     fprintf(stdout,"max corr = %f\n\n", maxval);
00528 
00529 
00530 
00531     /* note time at beginning of actual registration process */
00532     gettimeofday(&findgain, NULL);
00533 
00534 
00535 
00536     FINDGAINADJ(img1, img2, *x_off, *y_off, gain_adj, bias_adj);
00537         
00538 
00539     /* note finishing time */
00540     gettimeofday(&end, NULL);
00541 
00542 
00543 
00544     the_end:
00545 
00546     if (NAMES) {
00547             MESSAGE('I',"");
00548         sprintf(msg," Computed x-offset: %d",*x_off);
00549             MESSAGE('I',msg);
00550         sprintf(msg," Computed y-offset: %d",*y_off);
00551             MESSAGE('I',msg);
00552             MESSAGE('I',"");
00553         sprintf(msg,"Time to compute first correlation = %ld ms.", delay(find_first_reg, findreg));
00554             MESSAGE('I',msg);
00555         sprintf(msg,"Time to compute precise correlation = %ld ms.", delay(findreg, findgain));
00556             MESSAGE('I',msg);
00557         sprintf(msg,"Time to compute gain/bias registration = %ld ms.", delay(findgain, end));
00558             MESSAGE('I',msg);
00559         sprintf(msg,"Total time required = %ld ms.", delay(start, end));
00560             MESSAGE('I',msg);
00561         MESSAGE('I'," ...............");
00562     }
00563 
00564     /* if (TIMES) TIMING(T_EXIT); */
00565 
00566 }