IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 27672 for trunk/psModules


Ignore:
Timestamp:
Apr 13, 2010, 4:42:21 PM (16 years ago)
Author:
eugene
Message:

modify pmFootprintCullPeaks to reduce repeated allocs

Location:
trunk/psModules/src
Files:
2 added
7 edited

Legend:

Unmodified
Added
Removed
  • trunk/psModules/src/objects/Makefile.am

    r27657 r27672  
    1212        pmFootprintCullPeaks.c \
    1313        pmFootprintFind.c \
     14        pmFootprintSpans.c \
    1415        pmFootprintFindAtPoint.c \
    1516        pmFootprintIDs.c \
     
    7576        pmSpan.h \
    7677        pmFootprint.h \
     78        pmFootprintSpans.h \
    7779        pmPeaks.h \
    7880        pmMoments.h \
  • trunk/psModules/src/objects/pmFootprint.c

    r23187 r27672  
    5353    assert(nspan >= 0);
    5454    footprint->npix = 0;
     55    footprint->nspans = 0; // we may allocate more spans than we set -- this is the number of active spans
    5556    footprint->spans = psArrayAllocEmpty(nspan);
    5657    footprint->peaks = psArrayAlloc(0);
     
    7374}
    7475
     76bool pmFootprintAllocEmptySpans (pmFootprint *footprint, int nSpans) {
     77
     78    psArrayRealloc (footprint->spans, nSpans);
     79    for (int i = 0; i < nSpans; i++) {
     80        footprint->spans->data[i] = pmSpanAlloc(-1, -1, -1);
     81    }
     82    footprint->spans->n = nSpans;
     83    return true;
     84}
     85
     86// reset the footprint containers
     87bool pmFootprintInit(pmFootprint *footprint) {
     88
     89    footprint->bbox.x0 = footprint->bbox.y0 = 0;
     90    footprint->bbox.x1 = footprint->bbox.y1 = -1;
     91    footprint->nspans = 0;
     92    return true;
     93}
     94
    7595bool pmFootprintTest(const psPtr ptr) {
    7696    return (psMemGetDeallocator(ptr) == (psFreeFunc)footprintFree);
     
    103123    psArrayAdd(fp->spans, 1, sp);
    104124    psFree(sp);
    105    
     125
     126    fp->nspans ++;
     127
    106128    fp->npix += x1 - x0 + 1;
    107129
    108     if (fp->spans->n == 1) {
     130    if (fp->nspans == 1) {
    109131        fp->bbox.x0 = x0;
    110132        fp->bbox.x1 = x1;
     
    121143}
    122144
     145
     146// Set the next available elements of the nSpan entry in footprint->spans
     147pmSpan *pmFootprintSetSpan(pmFootprint *fp,     // the footprint to add to
     148                           const int y,         // row to add
     149                           int x0,              // range of
     150                           int x1) {            // columns
     151
     152    if (x1 < x0) {
     153        int tmp = x0;
     154        x0 = x1;
     155        x1 = tmp;
     156    }
     157
     158    int N = fp->nspans;
     159    if (N == fp->spans->n) {
     160        // if we need more space, extend fp->spans as needed
     161        int Nalloc = fp->spans->n + 100;
     162        psArrayRealloc(fp->spans, Nalloc);
     163        fp->spans->n = Nalloc;
     164        for (int i = N; i < fp->spans->n; i++) {
     165            fp->spans->data[i] = pmSpanAlloc(-1, -1, -1);
     166        }
     167    }
     168
     169    pmSpan *span = fp->spans->data[N];
     170    span->y = y;
     171    span->x0 = x0;
     172    span->x1 = x1;
     173
     174    fp->nspans ++;
     175
     176    fp->npix += x1 - x0 + 1;
     177
     178    if (fp->nspans == 1) {
     179        fp->bbox.x0 = x0;
     180        fp->bbox.x1 = x1;
     181        fp->bbox.y0 = y;
     182        fp->bbox.y1 = y;
     183    } else {
     184        if (x0 < fp->bbox.x0) fp->bbox.x0 = x0;
     185        if (x1 > fp->bbox.x1) fp->bbox.x1 = x1;
     186        if (y < fp->bbox.y0) fp->bbox.y0 = y;
     187        if (y > fp->bbox.y1) fp->bbox.y1 = y;
     188    }
     189
     190    return span;
     191}
     192
    123193void pmFootprintSetBBox(pmFootprint *fp) {
    124194    assert (fp != NULL);
    125     if (fp->spans->n == 0) {
     195    if (fp->nspans == 0) {
    126196        return;
    127197    }
     
    132202    int y1 = sp->y;
    133203
    134     for (int i = 1; i < fp->spans->n; i++) {
     204    for (int i = 1; i < fp->nspans; i++) {
    135205        sp = fp->spans->data[i];
    136206       
     
    150220   assert (fp != NULL);
    151221   int npix = 0;
    152    for (int i = 0; i < fp->spans->n; i++) {
     222   for (int i = 0; i < fp->nspans; i++) {
    153223       pmSpan *span = fp->spans->data[i];
    154224       npix += span->x1 - span->x0 + 1;
  • trunk/psModules/src/objects/pmFootprint.h

    r24875 r27672  
    1313#include <pslib.h>
    1414#include "pmSpan.h"
    15 
     15#include "pmFootprintSpans.h"
    1616
    1717typedef struct {
    1818    const int id;                       //!< unique ID
    1919    int npix;                           //!< number of pixels in this pmFootprint
    20     psArray *spans;                     //!< the pmSpans
     20    int nspans;
     21    psArray *spans;                     //!< the allocated pmSpans
    2122    psRegion bbox;                      //!< the pmFootprint's bounding box
    2223    psArray *peaks;                     //!< the peaks lying in this footprint
     
    2627
    2728pmFootprint *pmFootprintAlloc(int nspan, const psImage *img);
     29bool pmFootprintInit(pmFootprint *footprint);
    2830bool pmFootprintTest(const psPtr ptr);
     31
     32bool pmFootprintAllocEmptySpans (pmFootprint *footprint, int nSpans);
    2933
    3034pmFootprint *pmFootprintNormalize(pmFootprint *fp);
     
    3741                           int x1);    //          columns
    3842
     43pmSpan *pmFootprintSetSpan(pmFootprint *fp,     // the footprint to add to
     44                           const int y,         // row to add
     45                           int x0,              // range of
     46                           int x1);             // columns
     47
    3948psArray *pmFootprintsFind(const psImage *img, const float threshold, const int npixMin);
    40 pmFootprint *pmFootprintsFindAtPoint(const psImage *img,
    41                                     const float threshold,
    42                                     const psArray *peaks,
    43                                     int row, int col);
     49
     50bool pmFootprintsFindAtPoint(pmFootprint *fp,
     51                             pmFootprintSpans *fpSpans,
     52                             const psImage *img,     // image to search
     53                             psImage *mask,
     54                             const float threshold,   // Threshold
     55                             const psArray *peaks, // array of peaks; finding one terminates search for footprint
     56                             int row, int col);
     57
     58// pmFootprint *pmFootprintsFindAtPoint(const psImage *img,
     59//                                     const float threshold,
     60//                                     const psArray *peaks,
     61//                                     int row, int col);
     62
     63bool pmFootprintSpansBuild(pmFootprint *fp, // the footprint that we're building
     64                           pmFootprintSpans *fpSpans,
     65                           const psImage *img, // the psImage we're working on
     66                           psImage *mask, // the associated masks
     67                           const float threshold // Threshold
     68    );
    4469
    4570psArray *pmFootprintArrayGrow(const psArray *footprints, int r);
  • trunk/psModules/src/objects/pmFootprintCullPeaks.c

    r26893 r27672  
    2020#include "pmSpan.h"
    2121#include "pmFootprint.h"
     22#include "pmFootprintSpans.h"
    2223#include "pmPeaks.h"
     24
     25bool dumpfootprints (pmFootprint *fp, pmFootprintSpans *fpSp);
    2326
    2427 /*
     
    6568    psArrayAdd (brightPeaks, 128, fp->peaks->data[0]);
    6669
     70    // allocate the peakFootprint and peakFPSpans containers -- these are re-used by pmFootprintsFindAtPoint to minimize allocs in this function
     71    pmFootprint *peakFootprint = pmFootprintAlloc(fp->nspans, subImg);
     72    pmFootprintSpans *peakFPSpans = pmFootprintSpansAlloc(2*fp->nspans);
     73
     74    // allocate empty spans for the footprints
     75    pmFootprintAllocEmptySpans(peakFootprint, fp->nspans);
     76
     77    psImage *subMask = psImageCopy(NULL, subImg, PS_TYPE_IMAGE_MASK);
     78
    6779    // The brightest peak is always safe; go through other peaks trying to cull them
    6880    for (int i = 1; i < fp->peaks->n; i++) { // n.b. fp->peaks->n can change within the loop
     
    98110        }
    99111
    100         // XXX this is a bit expensive: psImageAlloc for every peak contained in this footprint
    101         // perhaps this should alloc a single ID image above and pass it in to be set.
    102 
    103         // XXX optionally use the faster pmFootprintsFind if the subimage size is large (eg, M31)
    104 
    105         // XXX this is not quite there yet:
    106         // pmFootprint *peakFootprint = NULL;
    107         // int area = subImg->numCols * subImg->numRows;
    108         // if (area > 30000) {
    109 
    110         pmFootprint *peakFootprint = pmFootprintsFindAtPoint(subImg, threshold, brightPeaks, peak->y, peak->x);
     112        // init peakFootprint here?
     113        pmFootprintsFindAtPoint(peakFootprint, peakFPSpans, subImg, subMask, threshold, brightPeaks, peak->y, peak->x);
     114        if (peakFPSpans->nStartSpans > 2000) {
     115            // dumpfootprints(peakFootprint, peakFPSpans);
     116            // fprintf (stderr, "big footprint %d : %d\n", peakFootprint->nspans, peakFPSpans->nStartSpans);
     117            fprintf (stderr, "big test footprint: %f %f to %f %f (%d pix)\n", peakFootprint->bbox.x0, peakFootprint->bbox.y0, peakFootprint->bbox.x1, peakFootprint->bbox.y1, peakFootprint->npix);
     118        }
    111119
    112120        // at this point brightPeaks only has the peaks brighter than the current
    113121
    114         // XXX need to supply the image here
    115122        // we set the IDs to either 1 (in peak) or 0 (not in peak)
    116123        pmSetFootprintID (idImg, peakFootprint, IN_PEAK);
    117         psFree(peakFootprint);
    118124
    119125        // If this peak has not already been assigned to a source, then we can look for any
     
    144150    psFree(subImg);
    145151    psFree(subWt);
     152    psFree(subMask);
     153    psFree(peakFootprint);
     154    psFree(peakFPSpans);
    146155
    147156    return PS_ERR_NONE;
     
    149158
    150159
    151  /*
    152   * Examine the peaks in a pmFootprint, and throw away the ones that are not sufficiently
    153   * isolated.  More precisely, for each peak find the highest coll that you'd have to traverse
    154   * to reach a still higher peak --- and if that coll's more than nsigma DN below your
    155   * starting point, discard the peak.
    156   */
    157 psErrorCode pmFootprintCullPeaks_OLD(const psImage *img, // the image wherein lives the footprint
    158                                  const psImage *weight, // corresponding variance image
    159                                  pmFootprint *fp, // Footprint containing mortal peaks
    160                                  const float nsigma_delta, // how many sigma above local background a peak needs to be to survive
    161                                  const float fPadding, // fractional padding added to stdev since bright peaks have unreasonably high significance
    162                                  const float min_threshold) { // minimum permitted coll height
    163     assert (img != NULL); assert (img->type.type == PS_TYPE_F32);
    164     assert (weight != NULL); assert (weight->type.type == PS_TYPE_F32);
    165     assert (img->row0 == weight->row0 && img->col0 == weight->col0);
    166     assert (fp != NULL);
     160bool dumpfootprints (pmFootprint *fp, pmFootprintSpans *fpSp) {
    167161
    168     if (fp->peaks == NULL || fp->peaks->n == 0) { // nothing to do
    169         return PS_ERR_NONE;
     162    FILE *f1 = fopen ("fp.dat", "w");
     163    if (!f1) return false;
     164
     165    for (int i = 0; i < fp->nspans; i++) {
     166        pmSpan *sp = fp->spans->data[i];
     167        fprintf (f1, "%d - %d : %d\n", sp->x0, sp->x1, sp->y);
    170168    }
     169    fclose (f1);
    171170
    172     psRegion subRegion;                 // desired subregion; 1 larger than bounding box (grr)
    173     subRegion.x0 = fp->bbox.x0; subRegion.x1 = fp->bbox.x1 + 1;
    174     subRegion.y0 = fp->bbox.y0; subRegion.y1 = fp->bbox.y1 + 1;
    175     const psImage *subImg = psImageSubset((psImage *)img, subRegion);
    176     const psImage *subWt = psImageSubset((psImage *)weight, subRegion);
    177     assert (subImg != NULL && subWt != NULL);
    178     //
    179     // We need a psArray of peaks brighter than the current peak.  We'll fake this
    180     // by reusing the fp->peaks but lying about n.
    181     //
    182     // We do this for efficiency (otherwise I'd need two peaks lists), and we are
    183     // rather too chummy with psArray in consequence.  But it works.
    184     //
    185     psArray *brightPeaks = psArrayAlloc(0);
    186     psFree(brightPeaks->data);
    187     brightPeaks->data = psMemIncrRefCounter(fp->peaks->data);// use the data from fp->peaks
    188     //
    189     // The brightest peak is always safe; go through other peaks trying to cull them
    190     //
    191     for (int i = 1; i < fp->peaks->n; i++) { // n.b. fp->peaks->n can change within the loop
    192         const pmPeak *peak = fp->peaks->data[i];
    193         int x = peak->x - subImg->col0;
    194         int y = peak->y - subImg->row0;
    195         //
    196         // Find the level nsigma below the peak that must separate the peak
    197         // from any of its friends
    198         //
    199         assert (x >= 0 && x < subImg->numCols && y >= 0 && y < subImg->numRows);
     171    FILE *f2 = fopen ("fpSp.dat", "w");
     172    if (!f2) return false;
    200173
    201         // stdev ~ sqrt(flux) : for bright regions (f ~ 1e6, sqrt(f) = 1e3 -- unrealistically small deviations)
    202         // add additional padding in quadrature:
    203 
    204         const float stdev = sqrt(subWt->data.F32[y][x]);
    205         const float flux = subImg->data.F32[y][x];
    206         const float fStdev = stdev/flux;
    207         const float stdev_pad = flux * hypot(fStdev, fPadding);
    208         float threshold = flux - nsigma_delta*stdev_pad;
    209 
    210         if (isnan(threshold) || threshold < min_threshold) {
    211 #if 1     // min_threshold is assumed to be below the detection threshold,
    212           // so all the peaks are pmFootprint, and this isn't the brightest
    213             // XXX mark peak to be dropped
    214             (void)psArrayRemoveIndex(fp->peaks, i);
    215             i--;                        // we moved everything down one
    216             continue;
    217 #else
    218 #error n.b. We will be running LOTS of checks at this threshold, so only find the footprint once
    219             threshold = min_threshold;
    220 #endif
    221         }
    222 
    223         // XXX EAM : if stdev >= 0, i'm not sure how this can ever be true?
    224         if (threshold > subImg->data.F32[y][x]) {
    225             threshold = subImg->data.F32[y][x] - 10*FLT_EPSILON;
    226         }
    227 
    228         // XXX this is a bit expensive: psImageAlloc for every peak contained in this footprint
    229         // perhaps this should alloc a single ID image above and pass it in to be set.
    230 
    231         const int peak_id = 1;          // the ID for the peak of interest
    232         brightPeaks->n = i;             // only stop at a peak brighter than we are
    233 
    234         // XXX optionally use the faster pmFootprintsFind if the subimage size is large (eg, M31)
    235 
    236         pmFootprint *peakFootprint = pmFootprintsFindAtPoint(subImg, threshold, brightPeaks, peak->y, peak->x);
    237         brightPeaks->n = 0;             // don't double free
    238         psImage *idImg = pmSetFootprintID(NULL, peakFootprint, peak_id);
    239         psFree(peakFootprint);
    240 
    241         // Check if any of the previous (brighter) peaks are within the footprint of this peak
    242         // If so, the current peak is bogus; drop it.
    243         int j;
    244         for (j = 0; j < i; j++) {
    245             const pmPeak *peak2 = fp->peaks->data[j];
    246             int x2 = peak2->x - subImg->col0;
    247             int y2 = peak2->y - subImg->row0;
    248             const int peak2_id = idImg->data.S32[y2][x2]; // the ID for some other peak
    249 
    250             if (peak2_id == peak_id) {  // There's a brighter peak within the footprint above
    251                 ;                       // threshold; so cull our initial peak
    252                 (void)psArrayRemoveIndex(fp->peaks, i);
    253                 i--;                    // we moved everything down one
    254                 break;
    255             }
    256         }
    257         if (j == i) {
    258             j++;
    259         }
    260 
    261         psFree(idImg);
     174    for (int i = 0; i < fpSp->nStartSpans; i++) {
     175        pmStartSpan *sp = fpSp->startspans->data[i];
     176        if (!sp->span) continue;
     177        fprintf (f2, "%d - %d : %d\n", sp->span->x0, sp->span->x1, sp->span->y);
    262178    }
    263 
    264     brightPeaks->n = 0; psFree(brightPeaks);
    265     psFree(subImg);
    266     psFree(subWt);
    267 
    268     return PS_ERR_NONE;
     179    fclose (f2);
     180    return true;
    269181}
    270 
  • trunk/psModules/src/objects/pmFootprintFindAtPoint.c

    r26893 r27672  
    2121#include "pmFootprint.h"
    2222#include "pmPeaks.h"
    23 
    24 /*
    25  * A data structure to hold the starting point for a search for pixels above threshold,
    26  * used by pmFootprintsFindAtPoint
    27  *
    28  * We don't want to find this span again --- it's already part of the footprint ---
    29  * so we set appropriate mask bits
    30  *
    31  * EAM : these function were confusingly using "startspan" and "spartspan"
    32  * I've rationalized them all to 'startspan'
    33  */
    34 
    35 //
    36 // An enum for what we should do with a Startspan
    37 //
    38 typedef enum {PM_SSPAN_DOWN = 0,        // scan down from this span
    39               PM_SSPAN_UP,              // scan up from this span
    40               PM_SSPAN_RESTART,         // restart scanning from this span
    41               PM_SSPAN_DONE             // this span is processed
    42 } PM_SSPAN_DIR;                         // How to continue searching
    43 //
    44 // An enum for mask's pixel values.  We're looking for pixels that are above threshold, and
    45 // we keep extra book-keeping information in the PM_SSPAN_STOP plane.  It's simpler to be
    46 // able to check for
    47 //
    48 enum {
    49     PM_SSPAN_INITIAL = 0x0,             // initial state of pixels.
    50     PM_SSPAN_DETECTED = 0x1,            // we've seen this pixel
    51     PM_SSPAN_STOP = 0x2                 // you may stop searching when you see this pixel
    52 };
    53 //
    54 // The struct that remembers how to [re-]start scanning the image for pixels
    55 //
    56 typedef struct {
    57     const pmSpan *span;                 // save the pixel range
    58     PM_SSPAN_DIR direction;             // How to continue searching
    59     bool stop;                          // should we stop searching?
    60 } Startspan;
    61 
    62 static void startspanFree(Startspan *sspan) {
    63     psFree(sspan->span);
    64 }
    65 
    66 static Startspan *
    67 StartspanAlloc(const pmSpan *span,      // The span in question
    68                psImage *mask,           // Pixels that we've already detected
    69                const PM_SSPAN_DIR dir   // Should we continue searching towards the top of the image?
    70     ) {
    71     Startspan *sspan = psAlloc(sizeof(Startspan));
    72     psMemSetDeallocator(sspan, (psFreeFunc)startspanFree);
    73 
    74     sspan->span = psMemIncrRefCounter((void *)span);
    75     sspan->direction = dir;
    76     sspan->stop = false;
    77 
    78     if (mask != NULL) {                 // remember that we've detected these pixels
    79         psImageMaskType *mpix = &mask->data.PS_TYPE_IMAGE_MASK_DATA[span->y - mask->row0][span->x0 - mask->col0];
    80 
    81         for (int i = 0; i <= span->x1 - span->x0; i++) {
    82             mpix[i] |= PM_SSPAN_DETECTED;
    83             if (mpix[i] & PM_SSPAN_STOP) {
    84                 sspan->stop = true;
    85             }
    86         }
    87     }
    88 
    89     return sspan;
    90 }
    91 
    92 //
    93 // Add a new Startspan to an array of Startspans.  Iff we see a stop bit, return true
    94 //
    95 static bool add_startspan(psArray *startspans, // the saved Startspans
    96                           const pmSpan *sp, // the span in question
    97                           psImage *mask, // mask of detected/stop pixels
    98                           const PM_SSPAN_DIR dir) { // the desired direction to search
    99     if (dir == PM_SSPAN_RESTART) {
    100         if (add_startspan(startspans, sp, mask,  PM_SSPAN_UP) ||
    101             add_startspan(startspans, sp, NULL, PM_SSPAN_DOWN)) {
    102             return true;
    103         }
    104     } else {
    105         Startspan *sspan = StartspanAlloc(sp, mask, dir);
    106         if (sspan->stop) {              // we detected a stop bit
    107             psFree(sspan);              // don't allocate new span
    108 
    109             return true;
    110         } else {
    111             psArrayAdd(startspans, 1, sspan);
    112             psFree(sspan);              // as it's now owned by startspans
    113         }
    114     }
    115 
    116     return false;
    117 }
     23#include "pmFootprintSpans.h"
    11824
    11925/************************************************************************************************************/
    12026/*
    121  * Search the image for pixels above threshold, starting at a single Startspan.
     27 * Search the image for pixels above threshold, starting at a single pmStartSpan.
    12228 * We search the array looking for one to process; it'd be better to move the
    12329 * ones that we're done with to the end, but it probably isn't worth it for
     
    12632 * This is the guts of pmFootprintsFindAtPoint
    12733 */
    128 static bool do_startspan(pmFootprint *fp, // the footprint that we're building
    129                          const psImage *img, // the psImage we're working on
    130                          psImage *mask, // the associated masks
    131                          const float threshold, // Threshold
    132                          psArray *startspans) { // specify which span to process next
     34bool pmFootprintSpansBuild(pmFootprint *fp, // the footprint that we're building
     35                           pmFootprintSpans *fpSpans,
     36                           const psImage *img, // the psImage we're working on
     37                           psImage *mask, // the associated masks
     38                           const float threshold // Threshold
     39    ) {
    13340    bool F32 = false;                   // is this an F32 image?
    13441    if (img->type.type == PS_TYPE_F32) {
     
    15259    /********************************************************************************************************/
    15360
    154     Startspan *sspan = NULL;
    155     for (int i = 0; i < startspans->n; i++) {
    156         sspan = startspans->data[i];
    157         if (sspan->direction != PM_SSPAN_DONE) {
     61    pmStartSpan *startspan = NULL;
     62    for (int i = 0; i < fpSpans->nStartSpans; i++) {
     63        startspan = fpSpans->startspans->data[i];
     64        if (startspan->direction != PM_STARTSPAN_DONE) {
    15865            break;
    15966        }
    160         if (sspan->stop) {
     67        if (startspan->stop) {
    16168            break;
    16269        }
    16370    }
    164     if (sspan == NULL || sspan->direction == PM_SSPAN_DONE) { // no more Startspans to process
     71    if (startspan == NULL || startspan->direction == PM_STARTSPAN_DONE) { // no more pmStartSpans to process
    16572        return false;
    16673    }
    167     if (sspan->stop) {                  // they don't want any more spans processed
     74    if (startspan->stop) {                  // they don't want any more spans processed
    16875        return false;
    16976    }
     77
    17078    /*
    17179     * Work
    17280     */
    173     const PM_SSPAN_DIR dir = sspan->direction;
     81    const PM_STARTSPAN_DIR dir = startspan->direction;
    17482    /*
    17583     * Set initial span to the startspan
    17684     */
    177     int x0 = sspan->span->x0 - col0, x1 = sspan->span->x1 - col0;
     85    int x0 = startspan->span->x0 - col0, x1 = startspan->span->x1 - col0;
    17886    /*
    17987     * Go through image identifying objects
    18088     */
    18189    int nx0, nx1 = -1;                  // new values of x0, x1
    182     const int di = (dir == PM_SSPAN_UP) ? 1 : -1; // how much i changes to get to the next row
     90    const int di = (dir == PM_STARTSPAN_UP) ? 1 : -1; // how much i changes to get to the next row
    18391    bool stop = false;                  // should I stop searching for spans?
    18492
    185     for (int i = sspan->span->y -row0 + di; i < numRows && i >= 0; i += di) {
     93    for (int i = startspan->span->y - row0 + di; i < numRows && i >= 0; i += di) {
    18694        imgRowF32 = img->data.F32[i];   // only one of
    18795        imgRowS32 = img->data.S32[i];   //      these is valid!
     
    195103        for (int j = x0 - 1; j >= -1; j--) {
    196104            double pixVal = (j < 0) ? threshold - 100 : (F32 ? imgRowF32[j] : imgRowS32[j]);
    197             if ((maskRow[j] & PM_SSPAN_DETECTED) || pixVal < threshold) {
     105            if ((maskRow[j] & PM_STARTSPAN_DETECTED) || pixVal < threshold) {
    198106                if (j < x0 - 1) {       // we found some pixels above threshold
    199107                    nx0 = j + 1;
     
    209117            // Search right in leftmost span
    210118            //
    211             //nx1 = 0;                  // make gcc happy
    212119            for (int j = nx0 + 1; j <= numCols; j++) {
    213120                double pixVal = (j >= numCols) ? threshold - 100 : (F32 ? imgRowF32[j] : imgRowS32[j]);
    214                 if ((maskRow[j] & PM_SSPAN_DETECTED) || pixVal < threshold) {
     121                if ((maskRow[j] & PM_STARTSPAN_DETECTED) || pixVal < threshold) {
    215122                    nx1 = j - 1;
    216123                    break;
     
    218125            }
    219126
    220             const pmSpan *sp = pmFootprintAddSpan(fp, i + row0, nx0 + col0, nx1 + col0);
    221 
    222             if (add_startspan(startspans, sp, mask, PM_SSPAN_RESTART)) {
     127            pmSpan *sp = pmFootprintSetSpan(fp, i + row0, nx0 + col0, nx1 + col0);
     128            bool status = pmFootprintSpansSet(fpSpans, sp, mask, PM_STARTSPAN_RESTART);
     129            // fprintf (stderr, "set 1: %d vs %d\n", fp->nspans, fpSpans->nStartSpans);
     130            if (status) {
    223131                stop = true;
    224132                break;
     
    238146        for (int j = nx1 + 1; j <= x1 + 1; j++) {
    239147            double pixVal = (j >= numCols) ? threshold - 100 : (F32 ? imgRowF32[j] : imgRowS32[j]);
    240             if (!(maskRow[j] & PM_SSPAN_DETECTED) && pixVal >= threshold) {
     148            if (!(maskRow[j] & PM_STARTSPAN_DETECTED) && pixVal >= threshold) {
    241149                int sx0 = j++;          // span that we're working on is sx0:sx1
    242150                int sx1 = -1;           // We know that if we got here, we'll also set sx1
    243151                for (; j <= numCols; j++) {
    244152                    double pixVal = (j >= numCols) ? threshold - 100 : (F32 ? imgRowF32[j] : imgRowS32[j]);
    245                     if ((maskRow[j] & PM_SSPAN_DETECTED) || pixVal < threshold) { // end of span
     153                    if ((maskRow[j] & PM_STARTSPAN_DETECTED) || pixVal < threshold) { // end of span
    246154                        sx1 = j;
    247155                        break;
     
    250158                assert (sx1 >= 0);
    251159
    252                 const pmSpan *sp;
     160                pmSpan *sp;
    253161                if (first) {
    254162                    if (sx1 <= x1) {
    255                         sp = pmFootprintAddSpan(fp, i + row0, sx0 + col0, sx1 + col0 - 1);
    256                         if (add_startspan(startspans, sp, mask, PM_SSPAN_DONE)) {
     163                        sp = pmFootprintSetSpan(fp, i + row0, sx0 + col0, sx1 + col0 - 1);
     164                        bool status = pmFootprintSpansSet(fpSpans, sp, mask, PM_STARTSPAN_DONE);
     165                        // fprintf (stderr, "set 2: %d vs %d\n", fp->nspans, fpSpans->nStartSpans);
     166                        if (status) {
    257167                            stop = true;
    258168                            break;
    259169                        }
    260170                    } else {            // overhangs to right
    261                         sp = pmFootprintAddSpan(fp, i + row0, sx0 + col0, x1 + col0);
    262                         if (add_startspan(startspans, sp, mask, PM_SSPAN_DONE)) {
     171                        sp = pmFootprintSetSpan(fp, i + row0, sx0 + col0, x1 + col0);
     172                        bool status = pmFootprintSpansSet(fpSpans, sp, mask, PM_STARTSPAN_DONE);
     173                        // fprintf (stderr, "set 3: %d vs %d\n", fp->nspans, fpSpans->nStartSpans);
     174                        if (status) {
    263175                            stop = true;
    264176                            break;
    265177                        }
    266                         sp = pmFootprintAddSpan(fp, i + row0, x1 + 1 + col0, sx1 + col0 - 1);
    267                         if (add_startspan(startspans, sp, mask, PM_SSPAN_RESTART)) {
     178                        sp = pmFootprintSetSpan(fp, i + row0, x1 + 1 + col0, sx1 + col0 - 1);
     179                        status = pmFootprintSpansSet(fpSpans, sp, mask, PM_STARTSPAN_RESTART);
     180                        // fprintf (stderr, "set 4: %d vs %d\n", fp->nspans, fpSpans->nStartSpans);
     181                        if (status) {
    268182                            stop = true;
    269183                            break;
     
    272186                    first = false;
    273187                } else {
    274                     sp = pmFootprintAddSpan(fp, i + row0, sx0 + col0, sx1 + col0 - 1);
    275                     if (add_startspan(startspans, sp, mask, PM_SSPAN_RESTART)) {
     188                    sp = pmFootprintSetSpan(fp, i + row0, sx0 + col0, sx1 + col0 - 1);
     189                    bool status = pmFootprintSpansSet(fpSpans, sp, mask, PM_STARTSPAN_RESTART);
     190                    // fprintf (stderr, "set 5: %d vs %d\n", fp->nspans, fpSpans->nStartSpans);
     191                    if (status) {
    276192                        stop = true;
    277193                        break;
     
    291207     */
    292208
    293     sspan->direction = PM_SSPAN_DONE;
     209    startspan->direction = PM_STARTSPAN_DONE;
    294210    return stop ? false : true;
    295211}
     
    307223 * are not sorted by increasing y, x0, x1.  If this matters to you, call
    308224 * pmFootprintNormalize()
     225 *
     226 * The calling function must supply a footprint allocated with a reasonable amount of space
     227 *
    309228 */
    310 pmFootprint *
    311 pmFootprintsFindAtPoint(const psImage *img,     // image to search
    312                        const float threshold,   // Threshold
    313                        const psArray *peaks, // array of peaks; finding one terminates search for footprint
    314                        int row, int col) { // starting position (in img's parent's coordinate system)
    315    assert(img != NULL);
    316 
    317    bool F32 = false;                    // is this an F32 image?
    318    if (img->type.type == PS_TYPE_F32) {
    319        F32 = true;
    320    } else if (img->type.type == PS_TYPE_S32) {
    321        F32 = false;
    322    } else {                             // N.b. You can't trivially add more cases here; F32 is just a bool
    323        psError(PS_ERR_UNKNOWN, true, "Unsupported psImage type: %d", img->type.type);
    324        return NULL;
    325    }
    326    psF32 *imgRowF32 = NULL;             // row pointer if F32
    327    psS32 *imgRowS32 = NULL;             //  "   "   "  "  !F32
    328 
    329    const int row0 = img->row0;
    330    const int col0 = img->col0;
    331    const int numRows = img->numRows;
    332    const int numCols = img->numCols;
    333 /*
    334  * Is point in image, and above threshold?
    335  */
    336    row -= row0; col -= col0;
    337    if (row < 0 || row >= numRows ||
    338        col < 0 || col >= numCols) {
     229
     230bool pmFootprintsFindAtPoint(pmFootprint *fp,
     231                             pmFootprintSpans *fpSpans,
     232                             const psImage *img,     // image to search
     233                             psImage *mask,
     234                             const float threshold,   // Threshold
     235                             const psArray *peaks, // array of peaks; finding one terminates search for footprint
     236                             int row, int col) { // starting position (in img's parent's coordinate system)
     237    psAssert(img, "image must be supplied");
     238    psAssert(fp, "footprint must be supplied");
     239    psAssert(fpSpans, "footprint spans must be supplied");
     240
     241    bool F32 = false;                    // is this an F32 image?
     242    if (img->type.type == PS_TYPE_F32) {
     243        F32 = true;
     244    } else if (img->type.type == PS_TYPE_S32) {
     245        F32 = false;
     246    } else {                             // N.b. You can't trivially add more cases here; F32 is just a bool
     247        psError(PS_ERR_UNKNOWN, true, "Unsupported psImage type: %d", img->type.type);
     248        return false;
     249    }
     250    psF32 *imgRowF32 = NULL;             // row pointer if F32
     251    psS32 *imgRowS32 = NULL;             //  "   "   "  "  !F32
     252
     253    const int row0 = img->row0;
     254    const int col0 = img->col0;
     255    const int numRows = img->numRows;
     256    const int numCols = img->numCols;
     257
     258    /*
     259     * Is point in image, and above threshold?
     260     */
     261    row -= row0; col -= col0;
     262    if (row < 0 || row >= numRows ||
     263        col < 0 || col >= numCols) {
    339264        psError(PS_ERR_BAD_PARAMETER_VALUE, true,
    340265                "row/col == (%d, %d) are out of bounds [%d--%d, %d--%d]",
    341266                row + row0, col + col0, row0, row0 + numRows - 1, col0, col0 + numCols - 1);
    342        return NULL;
    343    }
    344 
    345    double pixVal = F32 ? img->data.F32[row][col] : img->data.S32[row][col];
    346    if (pixVal < threshold) {
    347        return pmFootprintAlloc(0, img);
    348    }
    349 
    350    pmFootprint *fp = pmFootprintAlloc(1 + img->numRows/10, img);
    351 /*
    352  * We need a mask for two purposes; to indicate which pixels are already detected,
    353  * and to store the "stop" pixels --- those that, once reached, should stop us
    354  * looking for the rest of the pmFootprint.  These are generally set from peaks.
    355  */
    356    psImage *mask = psImageAlloc(numCols, numRows, PS_TYPE_IMAGE_MASK);
    357    P_PSIMAGE_SET_ROW0(mask, row0);
    358    P_PSIMAGE_SET_COL0(mask, col0);
    359    psImageInit(mask, PM_SSPAN_INITIAL);
    360    //
    361    // Set stop bits from peaks list
    362    //
    363    assert (peaks == NULL || peaks->n == 0 || psMemCheckPeak(peaks->data[0]));
    364    if (peaks != NULL) {
    365        for (int i = 0; i < peaks->n; i++) {
    366            pmPeak *peak = peaks->data[i];
    367            mask->data.PS_TYPE_IMAGE_MASK_DATA[peak->y - mask->row0][peak->x - mask->col0] |= PM_SSPAN_STOP;
    368        }
    369    }
    370 /*
    371  * Find starting span passing through (row, col)
    372  */
    373    psArray *startspans = psArrayAllocEmpty(1); // spans where we have to restart the search
    374 
    375    imgRowF32 = img->data.F32[row];      // only one of
    376    imgRowS32 = img->data.S32[row];      //      these is valid!
    377    psImageMaskType *maskRow = mask->data.PS_TYPE_IMAGE_MASK_DATA[row];
    378    {
    379        int i;
    380        for (i = col; i >= 0; i--) {
    381            pixVal = F32 ? imgRowF32[i] : imgRowS32[i];
    382            if ((maskRow[i] & PM_SSPAN_DETECTED) || pixVal < threshold) {
    383                break;
    384            }
    385        }
    386        int i0 = i;
    387        for (i = col; i < numCols; i++) {
    388            pixVal = F32 ? imgRowF32[i] : imgRowS32[i];
    389            if ((maskRow[i] & PM_SSPAN_DETECTED) || pixVal < threshold) {
    390                break;
    391            }
    392        }
    393        int i1 = i;
    394        const pmSpan *sp = pmFootprintAddSpan(fp, row + row0, i0 + col0 + 1, i1 + col0 - 1);
    395 
    396        (void)add_startspan(startspans, sp, mask, PM_SSPAN_RESTART);
    397    }
    398    /*
    399     * Now workout from those Startspans, searching for pixels above threshold
    400     */
    401    while (do_startspan(fp, img, mask, threshold, startspans)) continue;
    402    /*
    403     * Cleanup
    404     */
    405    psFree(mask);
    406    psFree(startspans);                  // restores the image pixel
    407 
    408    return fp;                           // pmFootprint really
     267        return false;
     268    }
     269
     270    double pixVal = F32 ? img->data.F32[row][col] : img->data.S32[row][col];
     271    if (pixVal < threshold) {
     272        return true;
     273    }
     274
     275    /*
     276     * We need a mask for two purposes; to indicate which pixels are already detected,
     277     * and to store the "stop" pixels --- those that, once reached, should stop us
     278     * looking for the rest of the pmFootprint.  These are generally set from peaks.
     279     */
     280
     281    pmFootprintInit(fp);
     282    pmFootprintSpansInit(fpSpans);
     283    psImageInit(mask, PM_STARTSPAN_INITIAL);
     284
     285    // fprintf (stderr, "init: %d vs %d\n", fp->nspans, fpSpans->nStartSpans);
     286
     287    //
     288    // Set stop bits from peaks list
     289    //
     290    assert (peaks == NULL || peaks->n == 0 || psMemCheckPeak(peaks->data[0]));
     291    if (peaks != NULL) {
     292        for (int i = 0; i < peaks->n; i++) {
     293            pmPeak *peak = peaks->data[i];
     294            mask->data.PS_TYPE_IMAGE_MASK_DATA[peak->y - mask->row0][peak->x - mask->col0] |= PM_STARTSPAN_STOP;
     295        }
     296    }
     297
     298    /*
     299     * Find starting span passing through (row, col)
     300     */
     301    imgRowF32 = img->data.F32[row];      // only one of
     302    imgRowS32 = img->data.S32[row];      //      these is valid!
     303    psImageMaskType *maskRow = mask->data.PS_TYPE_IMAGE_MASK_DATA[row];
     304    {
     305        int i;
     306        for (i = col; i >= 0; i--) {
     307            pixVal = F32 ? imgRowF32[i] : imgRowS32[i];
     308            if ((maskRow[i] & PM_STARTSPAN_DETECTED) || pixVal < threshold) {
     309                break;
     310            }
     311        }
     312        int i0 = i;
     313        for (i = col; i < numCols; i++) {
     314            pixVal = F32 ? imgRowF32[i] : imgRowS32[i];
     315            if ((maskRow[i] & PM_STARTSPAN_DETECTED) || pixVal < threshold) {
     316                break;
     317            }
     318        }
     319        int i1 = i;
     320        pmSpan *sp = pmFootprintSetSpan(fp, row + row0, i0 + col0 + 1, i1 + col0 - 1);
     321        pmFootprintSpansSet(fpSpans, sp, mask, PM_STARTSPAN_RESTART);
     322        // fprintf (stderr, "first: %d vs %d\n", fp->nspans, fpSpans->nStartSpans);
     323    }
     324    /*
     325     * Now workout from those pmStartSpans, searching for pixels above threshold
     326     */
     327    while (pmFootprintSpansBuild(fp, fpSpans, img, mask, threshold)) continue;
     328    /*
     329     * Cleanup
     330     */
     331    // psFree(mask);
     332    // psFree(startspans);                  // restores the image pixel
     333
     334    return fp;                           // pmFootprint really
    409335}
  • trunk/psModules/src/objects/pmFootprintIDs.c

    r20937 r27672  
    3232   const int row0 = fp->region.y0;
    3333
    34    for (int j = 0; j < fp->spans->n; j++) {
     34   for (int j = 0; j < fp->nspans; j++) {
    3535       const pmSpan *span = fp->spans->data[j];
    3636       psS32 *imgRow = idImage->data.S32[span->y - row0];
  • trunk/psModules/src/psmodules.h

    r27657 r27672  
    111111// the following headers are from psModule:objects
    112112#include <pmSpan.h>
     113#include <pmFootprintSpans.h>
    113114#include <pmFootprint.h>
    114115#include <pmPeaks.h>
Note: See TracChangeset for help on using the changeset viewer.