IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Nov 17, 2006, 1:00:31 PM (19 years ago)
Author:
magnier
Message:

significant cleanup, added bicube measure of fractional peak position

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/psModules/src/objects/pmPeaks.c

    r9881 r10048  
    66 *  @author EAM, IfA: significant modifications.
    77 *
    8  *  @version $Revision: 1.11 $ $Name: not supported by cvs2svn $
    9  *  @date $Date: 2006-11-07 09:07:42 $
     8 *  @version $Revision: 1.12 $ $Name: not supported by cvs2svn $
     9 *  @date $Date: 2006-11-17 23:00:31 $
    1010 *
    1111 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    2424
    2525/******************************************************************************
    26 myListAddPeak(): A private function which allocates a psArray, if the list
    27 argument is NULL, otherwise it adds the peak to that list.
    28 XXX EAM : changed the output to psArray
    29 XXX EAM : Switched row, col args
    30 XXX EAM : NOTE: this was changed in the call, so the new code is consistent
    31 *****************************************************************************/
    32 static psArray *myListAddPeak(psArray *list,
    33                               psS32 row,
    34                               psS32 col,
    35                               psF32 counts,
    36                               pmPeakType type)
    37 {
    38     psTrace("psModules.objects", 4, "---- %s() begin ----\n", __func__);
    39     pmPeak *tmpPeak = pmPeakAlloc(col, row, counts, type);
    40 
    41     if (list == NULL) {
    42         list = psArrayAllocEmpty(100);
    43     }
    44     psArrayAdd(list, 100, tmpPeak);
    45     psFree (tmpPeak);
    46     // XXX EAM : is this free appropriate?  (does psArrayAdd increment memory counter?)
    47 
    48     psTrace("psModules.objects", 4, "---- %s() end ----\n", __func__);
    49     return(list);
     26AddPeak(): A private function which allocates a psArray, if the peaks
     27argument is NULL, otherwise it adds the peak to that array.
     28XXX EAM : row,col now refer to image coords, NOT parent (since this is private)
     29XXX EAM : now also calculates fractional peak positions from 3x3 bicube region
     30*****************************************************************************/
     31static psArray *AddPeak(psArray *peaks,
     32                        const psImage *image,
     33                        psS32 row,
     34                        psS32 col,
     35                        pmPeakType type)
     36{
     37    psTrace(__func__, 5, "---- begin ----\n");
     38
     39    if (peaks == NULL) {
     40        peaks = psArrayAllocEmpty(100);
     41    }
     42
     43    // the peak position is in parent coordinates
     44    pmPeak *peak = pmPeakAlloc(col + image->col0, row + image->row0, image->data.F32[row][col], type);
     45
     46    // measure fractional peak position using the 3x3 bicube fit
     47
     48    // ix,iy must land on image with 1 pixel border
     49    int ix = PS_MAX (PS_MIN (col, image->numCols - 2), 1);
     50    int iy = PS_MAX (PS_MIN (row, image->numRows - 2), 1);
     51
     52    // calculate peak position relative to ix,iy
     53    psPolynomial2D *bicube = psImageBicubeFit (image, ix + image->col0, iy + image->row0);
     54    psPlane min = psImageBicubeMin (bicube);
     55    psFree (bicube);
     56
     57    // if min point is too deviant, use the peak value
     58    if ((fabs(min.x) < 1.5) && (fabs(min.y) < 1.5)) {
     59        peak->xf = min.x + ix + image->col0;
     60        peak->yf = min.y + iy + image->row0;
     61    } else {
     62        peak->xf = ix;
     63        peak->yf = iy;
     64    }
     65
     66    psArrayAdd(peaks, 100, peak);
     67    psFree (peak);
     68
     69    psTrace(__func__, 5, "---- end ----\n");
     70    return(peaks);
    5071}
    5172
    … …  
    306327    psArray *list = NULL;
    307328
    308     psU32 col0 = image->col0;
    309     psU32 row0 = image->row0;
    310 
    311     //
    312329    // Find peaks in row 0 only.
    313     //
    314330    row = 0;
    315331    tmpRow = getRowVectorFromImage((psImage *) image, row);
    … …  
    319335    for (psU32 i = 0 ; i < row1->n ; i++ ) {
    320336        col = row1->data.U32[i];
    321         //
    322         // Determine if pixel (0,0) is a peak.
    323         //
     337        // is pixel (0,0) is a peak?
    324338        if (col == 0) {
    325339            if ( (image->data.F32[row][col] >  image->data.F32[row][col+1]) &&
    … …  
    328342
    329343                if (image->data.F32[row][col] > threshold) {
    330                     list = myListAddPeak(list, row + row0, col + col0, image->data.F32[row][col], PM_PEAK_EDGE);
     344                    list = AddPeak(list, image, row, col, PM_PEAK_EDGE);
    331345                }
    332346            }
    … …  
    338352                    (image->data.F32[row][col] >= image->data.F32[row+1][col+1])) {
    339353                if (image->data.F32[row][col] > threshold) {
    340                     list = myListAddPeak(list, row + row0, col + col0, image->data.F32[row][col], PM_PEAK_EDGE);
     354                    list = AddPeak(list, image, row, col, PM_PEAK_EDGE);
    341355                }
    342356            }
    … …  
    347361                    (image->data.F32[row][col] >= image->data.F32[row+1][col-1])) {
    348362                if (image->data.F32[row][col] > threshold) {
    349                     list = myListAddPeak(list, row + row0, col + col0, image->data.F32[row][col], PM_PEAK_EDGE);
     363                    list = AddPeak(list, image, row, col, PM_PEAK_EDGE);
    350364                }
    351365            }
    … …  
    386400                        (image->data.F32[row][col] >= image->data.F32[row+1][col+1])) {
    387401                    myType = PM_PEAK_EDGE;
    388                     list = myListAddPeak(list, row + row0, col + col0, image->data.F32[row][col], myType);
     402                    list = AddPeak(list, image, row, col, myType);
    389403                }
    390404            } else if (col < (image->numCols - 1)) {
    … …  
    421435                        }
    422436
    423                         list = myListAddPeak(list, row + row0, col + col0, image->data.F32[row][col], myType);
    424 
    425                         # if (0)
    426 
    427                             psPolynomial2D *bicube = psImageBicubeFit (inSource->pixels, ix, iy);
    428                         psPlane min = psImageBicubeMin (bicube);
    429                         peak->xf = min.x;
    430                         peak->yf = min.y;
    431                         # endif
     437                        list = AddPeak(list, image, row, col, myType);
    432438
    433439                    }
    … …  
    442448                        (image->data.F32[row][col] >= image->data.F32[row+1][col])) {
    443449                    myType = PM_PEAK_EDGE;
    444                     list = myListAddPeak(list, row + row0, col + col0, image->data.F32[row][col], myType);
     450                    list = AddPeak(list, image, row, col, myType);
    445451                }
    446452            } else {
    … …  
    466472                    (image->data.F32[row][col] >  image->data.F32[row][col+1])) {
    467473                if (image->data.F32[row][col] > threshold) {
    468                     list = myListAddPeak(list, row + row0, col + col0, image->data.F32[row][col], PM_PEAK_EDGE);
     474                    list = AddPeak(list, image, row, col, PM_PEAK_EDGE);
    469475                }
    470476            }
    … …  
    476482                    (image->data.F32[row][col] >= image->data.F32[row][col+1])) {
    477483                if (image->data.F32[row][col] > threshold) {
    478                     list = myListAddPeak(list, row + row0, col + col0, image->data.F32[row][col], PM_PEAK_EDGE);
     484                    list = AddPeak(list, image, row, col, PM_PEAK_EDGE);
    479485                }
    480486            }
    … …  
    485491                    (image->data.F32[row][col] >  image->data.F32[row][col-1])) {
    486492                if (image->data.F32[row][col] > threshold) {
    487                     list = myListAddPeak(list, row + row0, col + col0, image->data.F32[row][col], PM_PEAK_EDGE);
     493                    list = AddPeak(list, image, row, col, PM_PEAK_EDGE);
    488494                }
    489495            }
Note: See TracChangeset for help on using the changeset viewer.