Index: /trunk/psphot/src/Makefile.am =================================================================== --- /trunk/psphot/src/Makefile.am (revision 13006) +++ /trunk/psphot/src/Makefile.am (revision 13007) @@ -18,4 +18,5 @@ libpsphot_la_SOURCES = \ psphotErrorCodes.c \ + pmFootprint.c \ psphotVersion.c \ psphotModelGroupInit.c \ @@ -62,4 +63,5 @@ include_HEADERS = \ psphot.h \ + pmFootprint.h \ psphotErrorCodes.h Index: /trunk/psphot/src/pmFootprint.c =================================================================== --- /trunk/psphot/src/pmFootprint.c (revision 13007) +++ /trunk/psphot/src/pmFootprint.c (revision 13007) @@ -0,0 +1,710 @@ +/*****************************************************************************/ +/* + * Handle pmFootprints + */ +#include +#include +#include +#include +#include "psphot.h" +#include "pmFootprint.h" + +static void spanFree(pmSpan *tmp) { + ; +} + +/* + * pmSpanAlloc() + */ +pmSpan *pmSpanAlloc(int y, int x0, int x1) +{ + pmSpan *span = (pmSpan *)psAlloc(sizeof(pmSpan)); + psMemSetDeallocator(span, (psFreeFunc) spanFree); + + span->y = y; + span->x0 = x0; + span->x1 = x1; + + return(span); +} + +bool pmIsSpan(const psPtr ptr) +{ + return (psMemGetDeallocator(ptr) == (psFreeFunc)spanFree); +} + +/************************************************************************************************************/ + +static void footprintFree(pmFootprint *tmp) +{ + if (!tmp) { + return; + } + + psTrace("psModules.objects", 5, "---- begin ----\n"); + + psFree(tmp->spans); + psFree(tmp->peaks); + + psTrace("psModules.objects", 5, "---- end ----\n"); +} + +/* + * pmFootprintAlloc(): Allocate the pmFootprint structure to NULL. + */ +pmFootprint *pmFootprintAlloc(int nspan, // number of spans expected in footprint + const psImage *image) // region footprint lives in +{ + psTrace("psModules.objects", 5, "---- begin ----\n"); + + static int id = 1; + pmFootprint *footprint = (pmFootprint *)psAlloc(sizeof(pmFootprint)); + *(int *)&footprint->id = id++; + psMemSetDeallocator(footprint, (psFreeFunc) footprintFree); + + assert(nspan >= 0); + footprint->npix = 0; + footprint->nspan = nspan; + footprint->spans = psArrayAllocEmpty(nspan); + footprint->peaks = psArrayAlloc(0); + + footprint->bbox.x0 = footprint->bbox.y0 = 0; + footprint->bbox.x1 = footprint->bbox.y1 = -1; + + if (image == NULL) { + footprint->region.x0 = footprint->region.y0 = 0; + footprint->region.x1 = footprint->region.y1 = -1; + } else { + footprint->region.x0 = image->col0; + footprint->region.x1 = image->col0 + image->numCols - 1; + footprint->region.y0 = image->row0; + footprint->region.y1 = image->row0 + image->numRows - 1; + } + + psTrace("psModules.objects", 5, "---- end ----\n"); + return(footprint); +} + +bool pmIsFootprint(const psPtr ptr) { + return (psMemGetDeallocator(ptr) == (psFreeFunc)footprintFree); +} + +// +// Add a span to a footprint, returning the number of pixels in the footprint +// +static int pmFootprintAddSpan(pmFootprint *fp, // the footprint to add to + const int y, // row to add + int x0, // range of + int x1) { // columns + + if (x1 < x0) { + int tmp = x0; + x0 = x1; + x1 = tmp; + } + + pmSpan *sp = pmSpanAlloc(y, x0, x1); + psArrayAdd(fp->spans, 1, sp); + psFree(sp); + + fp->npix += x1 - x0 + 1; + + if (fp->nspan == 1) { + fp->bbox.x0 = x0; + fp->bbox.x1 = x1; + fp->bbox.y0 = y; + fp->bbox.y1 = y; + } else { + if (x0 < fp->bbox.x0) fp->bbox.x0 = x0; + if (x1 > fp->bbox.x1) fp->bbox.x1 = x1; + if (y < fp->bbox.y0) fp->bbox.y0 = y; + if (y > fp->bbox.y1) fp->bbox.y1 = y; + } + + return fp->npix; +} + +void pmFootprintSetBBox(pmFootprint *fp) { + assert (fp != NULL); + if (fp->nspan == 0) { + return; + } + pmSpan *sp = fp->spans->data[0]; + int x0 = sp->x0; + int x1 = sp->x1; + int y0 = sp->y; + int y1 = sp->y; + + for (int i = 1; i < fp->spans->n; i++) { + sp = fp->spans->data[i]; + + if (sp->x0 < x0) x0 = sp->x0; + if (sp->x1 > x1) x1 = sp->x1; + if (sp->y < y0) y0 = sp->y; + if (sp->y > y1) y1 = sp->y; + } + + fp->bbox.x0 = x0; + fp->bbox.x1 = x1; + fp->bbox.y0 = y0; + fp->bbox.y1 = y1; +} + +int pmFootprintSetNpix(pmFootprint *fp) { + assert (fp != NULL); + int npix = 0; + for (int i = 0; i < fp->nspan; i++) { + pmSpan *span = fp->spans->data[i]; + npix += span->x1 - span->x0 + 1; + } + fp->npix = npix; + + return npix; +} + +/************************************************************************************************************/ + +typedef struct { /* run-length code for part of object*/ + int id; /* ID for object */ + int y; /* Row wherein WSPAN dwells */ + int x0, x1; /* inclusive range of columns */ +} WSPAN; + +/* + * comparison function for qsort; sort by ID then row + */ +static int +compar(const void *va, const void *vb) +{ + const WSPAN *a = va; + const WSPAN *b = vb; + + if(a->id < b->id) { + return(-1); + } else if(a->id > b->id) { + return(1); + } else { + return(a->y - b->y); + } +} + +/* + * Follow a chain of aliases, returning the final resolved value. + */ +static int +resolve_alias(const int *aliases, /* list of aliases */ + int id) /* alias to look up */ +{ + int resolved = id; /* resolved alias */ + + while(id != aliases[id]) { + resolved = id = aliases[id]; + } + + return(resolved); +} + +/* + * Go through an image, finding sets of connected pixels above threshold + * and assembling them into pmFootprints; the resulting set of objects + * is returned as a psArray + */ +psArray * +pmFindFootprints(const psImage *img, // image to search + const float threshold, // Threshold + const int npixMin) // minimum number of pixels in an acceptable pmFootprint +{ + int *aliases; /* aliases for object IDs */ + int *id_s; /* storage for id[cp] */ + int *idc, *idp; /* object IDs in current/previous row*/ + int i0; /* initial value of i */ + int id; /* object ID */ + int in_span; /* object ID of current WSPAN */ + int nspan = 0; /* number of spans */ + int nobj = 0; /* number of objects found */ + int x0 = 0; /* unpacked from a WSPAN */ + int size_aliases = 0; /* size of aliases[] array */ + int size_spans = 0; /* size of spans[] array */ + WSPAN *spans; /* row:x0,x1 for objects */ + int *tmp; /* used in swapping idc/idp */ + + assert(img != NULL); + + bool F32 = false; // is this an F32 image? + if (img->type.type == PS_TYPE_F32) { + F32 = true; + } else if (img->type.type == PS_TYPE_S32) { + F32 = false; + } else { // N.b. You can't trivially add more cases here; F32 is just a bool + psError(PS_ERR_UNKNOWN, true, "Unsupported psImage type: %d", img->type.type); + return NULL; + } + psF32 *imgRowF32 = NULL; // row pointer if F32 + psS32 *imgRowS32 = NULL; // " " " " !F32 + + const int row0 = img->row0; + const int col0 = img->col0; + const int numRows = img->numRows; + const int numCols = img->numCols; +/* + * Storage for arrays that identify objects by ID. We want to be able to + * refer to idp[-1] and idp[numCols], hence the (numCols + 2) + */ + id_s = psAlloc(2*(numCols + 2)*sizeof(int)); + memset(id_s, '\0', 2*(numCols + 2)*sizeof(int)); assert(id_s[0] == 0); + idc = id_s + 1; idp = idc + (numCols + 2); + + size_aliases = 1 + numRows/20; + aliases = psAlloc(size_aliases*sizeof(int)); + + size_spans = 1 + numRows/20; + spans = psAlloc(size_spans*sizeof(WSPAN)); +/* + * Go through image identifying objects + */ + for (int i = 0; i < numRows; i++) { + int j; + tmp = idc; idc = idp; idp = tmp; /* swap ID pointers */ + memset(idc, '\0', numCols*sizeof(int)); + + imgRowF32 = img->data.F32[i]; // only one of + imgRowS32 = img->data.S32[i]; // these is valid! + + in_span = 0; /* not in a span */ + for (j = 0; j < numCols; j++) { + double pixVal = F32 ? imgRowF32[j] : imgRowS32[j]; + if (pixVal < threshold) { + if (in_span) { + if(nspan >= size_spans) { + size_spans *= 2; + spans = psRealloc(spans, size_spans*sizeof(WSPAN)); + } + spans[nspan].id = in_span; + spans[nspan].y = i; + spans[nspan].x0 = x0; + spans[nspan].x1 = j - 1; + + nspan++; + + in_span = 0; + } + } else { /* a pixel to fix */ + if(idc[j - 1] != 0) { + id = idc[j - 1]; + } else if(idp[j - 1] != 0) { + id = idp[j - 1]; + } else if(idp[j] != 0) { + id = idp[j]; + } else if(idp[j + 1] != 0) { + id = idp[j + 1]; + } else { + id = ++nobj; + + if(id >= size_aliases) { + size_aliases *= 2; + aliases = psRealloc(aliases, size_aliases*sizeof(int)); + } + aliases[id] = id; + } + + idc[j] = id; + if(!in_span) { + x0 = j; in_span = id; + } +/* + * Do we need to merge ID numbers? If so, make suitable entries in aliases[] + */ + if(idp[j + 1] != 0 && idp[j + 1] != id) { + aliases[resolve_alias(aliases, idp[j + 1])] = + resolve_alias(aliases, id); + + idc[j] = id = idp[j + 1]; + } + } + } + + if(in_span) { + if(nspan >= size_spans) { + size_spans *= 2; + spans = psRealloc(spans, size_spans*sizeof(WSPAN)); + } + + assert(nspan < size_spans); /* we checked for space above */ + spans[nspan].id = in_span; + spans[nspan].y = i; + spans[nspan].x0 = x0; + spans[nspan].x1 = j - 1; + + nspan++; + } + } + + psFree(id_s); +/* + * Resolve aliases; first alias chains, then the IDs in the spans + */ + for (int i = 0; i < nspan; i++) { + spans[i].id = resolve_alias(aliases, spans[i].id); + } + + psFree(aliases); +/* + * Sort spans by ID, so we can sweep through them once + */ + if(nspan > 0) { + qsort(spans, nspan, sizeof(WSPAN), compar); + } +/* + * Build pmFootprints from the spans + */ + psArray *footprints = psArrayAlloc(nobj); + int n = 0; // number of pmFootprints + + if(nspan > 0) { + id = spans[0].id; + i0 = 0; + for (int i = 0; i <= nspan; i++) { /* nspan + 1 to catch last object */ + if(i == nspan || spans[i].id != id) { + pmFootprint *fp = pmFootprintAlloc(i - i0, img); + + for(; i0 < i; i0++) { + pmFootprintAddSpan(fp, spans[i0].y + row0, + spans[i0].x0 + col0, spans[i0].x1 + col0); + } + + if (fp->npix < npixMin) { + psFree(fp); + } else { + footprints->data[n++] = fp; + } + } + + id = spans[i].id; + } + } + + footprints = psArrayRealloc(footprints, n); +/* + * clean up + */ + psFree(spans); + + return footprints; +} + +/************************************************************************************************************/ +/* + * Worker routine for the pmSetFootprintArrayIds (and pmMergeFootprintArrays) + */ +static void +set_footprint_array_ids(psImage *idImage, + const psArray *footprints, // the footprints to insert + const bool relativeIDs) { // show IDs starting at 0, not pmFootprint->id + const pmFootprint *fp = footprints->data[0]; + const int col0 = fp->region.x0; + const int row0 = fp->region.y0; + + int id = 0; // first index will be 1 + for (int i = 0; i < footprints->n; i++) { + fp = footprints->data[i]; + if (relativeIDs) { + id++; + } else { + id = fp->id; + } + + for (int j = 0; j < fp->nspan; j++) { + const pmSpan *span = fp->spans->data[j]; + psS32 *imgRow = idImage->data.S32[span->y - row0]; + for(int k = span->x0 - col0; k <= span->x1 - col0; k++) { + imgRow[k] += id; + } + } + } +} + +/* + * Set an image to the value of footprint's ID whever they may fall + */ +psImage *pmSetFootprintArrayIDs(const psArray *footprints, // the footprints to insert + const bool relativeIDs) { // show IDs starting at 0, not pmFootprint->id + assert (footprints != NULL); + + if (footprints->n == 0) { + psError(PS_ERR_BAD_PARAMETER_SIZE, true, "You didn't provide any footprints"); + return NULL; + } + const pmFootprint *fp = footprints->data[0]; + assert(pmIsFootprint((const psPtr)fp)); + const int numCols = fp->region.x1 - fp->region.x0 + 1; + const int numRows = fp->region.y1 - fp->region.y0 + 1; + const int col0 = fp->region.x0; + const int row0 = fp->region.y0; + assert (numCols >= 0 && numRows >= 0); + + psImage *idImage = psImageAlloc(numCols, numRows, PS_TYPE_S32); + P_PSIMAGE_SET_ROW0(idImage, row0); + P_PSIMAGE_SET_COL0(idImage, col0); + psImageInit(idImage, 0); + /* + * do the work + */ + set_footprint_array_ids(idImage, footprints, relativeIDs); + + return idImage; + +} + +/************************************************************************************************************/ +/* + * Merge together two psArrays of pmFootprints neither of which is damaged. + * + * The returned psArray may contain elements of the inital psArrays (with + * their reference counters suitable incremented) + */ +psArray *pmMergeFootprintArrays(const psArray *footprints1, + const psArray *footprints2) { + assert (footprints1->n == 0 || pmIsFootprint(footprints1->data[0])); + assert (footprints2->n == 0 || pmIsFootprint(footprints2->data[0])); + + if (footprints1->n == 0 || footprints2->n == 0) { // nothing to do but put copies on merged + const psArray *old = (footprints1->n == 0) ? footprints2 : footprints1; + + psArray *merged = psArrayAllocEmpty(old->n); + for (int i = 0; i < old->n; i++) { + psArrayAdd(merged, 1, old->data[i]); + } + + return merged; + } + /* + * We have real work to do as some pmFootprints in footprints2 may overlap + * with footprints1 + */ + { + pmFootprint *fp1 = footprints1->data[0]; + pmFootprint *fp2 = footprints2->data[0]; + if (fp1->region.x0 != fp2->region.x0 || + fp1->region.x1 != fp2->region.x1 || + fp1->region.y0 != fp2->region.y0 || + fp1->region.y1 != fp2->region.y1) { + psError(PS_ERR_BAD_PARAMETER_SIZE, true, + "The two pmFootprint arrays correspnond to different-sized regions"); + return NULL; + } + } + /* + * We'll insert first one set of footprints then the other into an image, then + * extract a footprint from the result --- this is magically what we want. + */ + psImage *idImage = pmSetFootprintArrayIDs(footprints1, true); + set_footprint_array_ids(idImage, footprints2, true); + + psArray *merged = pmFindFootprints(idImage, 0.5, 1); + assert (merged != NULL); + psFree(idImage); + /* + * Now assign the peaks appropriately. We could do this more efficiently + * using idImage (which we just freed), but this is easy any probably fast enough + */ + const psArray *peaks; + peaks = pmFootprintArrayToPeaks(footprints1); + pmPeaksAssignToFootprints(merged, peaks); + peaks = pmFootprintArrayToPeaks(footprints2); + pmPeaksAssignToFootprints(merged, peaks); + + return merged; +} + +/************************************************************************************************************/ +/* + * Given a psArray of pmFootprints and another of pmPeaks, assign the peaks to the + * footprints in which that fall; if they _don't_ fall in a footprint, add a suitable + * one to the list. + */ +psErrorCode +pmPeaksAssignToFootprints(psArray *footprints, // the pmFootprints + const psArray *peaks) { // the pmPeaks + assert (footprints != NULL); + assert (footprints->n == 0 || pmIsFootprint(footprints->data[0])); + assert (peaks != NULL); + assert (peaks->n == 0 || pmIsPeak(peaks->data[0])); + + if (footprints->n == 0) { + if (peaks->n > 0) { + return psError(PS_ERR_BAD_PARAMETER_SIZE, true, "Your list of footprints is empty"); + } + return PS_ERR_NONE; + } + /* + * Create an image filled with the object IDs, and use it to assign pmPeaks to the + * objects + */ + psImage *ids = pmSetFootprintArrayIDs(footprints, true); + assert (ids != NULL); + assert (ids->type.type == PS_TYPE_S32); + const int row0 = ids->row0; + const int col0 = ids->col0; + const int numRows = ids->numRows; + const int numCols = ids->numCols; + + for (int i = 0; i < peaks->n; i++) { + pmPeak *peak = peaks->data[i]; + const int x = peak->x - col0; + const int y = peak->y - row0; + + assert (x >= 0 && x < numCols && y >= 0 && y < numRows); + int id = ids->data.S32[y][x - col0]; + + if (id == 0) { // peak isn't in a footprint, so make one for it + pmFootprint *nfp = pmFootprintAlloc(1, ids); + pmFootprintAddSpan(nfp, y, x, x); + psArrayAdd(footprints, 1, nfp); + id = footprints->n; + } + + assert (id >= 1 && id <= footprints->n); + pmFootprint *fp = footprints->data[id - 1]; + psArrayAdd(fp->peaks, 5, peak); + } + + psFree(ids); + // + // Make sure that peaks within each footprint are sorted + // + for (int i = 0; i < footprints->n; i++) { + pmFootprint *fp = footprints->data[i]; + fp->peaks = psArraySort (fp->peaks, pmPeakSortBySN); + } + + return PS_ERR_NONE; +} + +/************************************************************************************************************/ + /* + * Examine the peaks in a pmFootprint, and throw away the ones that are not sufficiently + * isolated. More precisely, for each peak find the highest coll that you'd have to traverse + * to reach a still higher peak --- and if that coll's more than nsigma DN below your + * starting point, discard the peak. + */ +psErrorCode pmFootprintCullPeaks(const psImage *img, // the image wherein lives the footprint + const psImage *weight, // corresponding variance image + pmFootprint *fp, // Footprint containing mortal peaks + const float nsigma_delta, // how many sigma above local background a peak + // needs to be to survive + const float min_threshold) { // minimum permitted coll height + assert (img != NULL); assert (img->type.type == PS_TYPE_F32); + assert (weight != NULL); assert (weight->type.type == PS_TYPE_F32); + assert (img->row0 == weight->row0 && img->col0 == weight->col0); + assert (fp != NULL); + + if (fp->peaks == NULL || fp->peaks->n == 0) { // nothing to do + return PS_ERR_NONE; + } + + psRegion subRegion; // desired subregion; 1 larger than bounding box (grr) + subRegion.x0 = fp->bbox.x0; subRegion.x1 = fp->bbox.x1 + 1; + subRegion.y0 = fp->bbox.y0; subRegion.y1 = fp->bbox.y1 + 1; + const psImage *subImg = psImageSubset((psImage *)img, subRegion); + const psImage *subWt = psImageSubset((psImage *)weight, subRegion); + assert (subImg != NULL && subWt != NULL); + // + // The brightest peak is always safe; go through other peaks trying to cull them + // + for (int i = 1; i < fp->peaks->n; i++) { // n.b. fp->peaks->n can change within the loop + const pmPeak *peak = fp->peaks->data[i]; + int x = peak->x - subImg->col0; + int y = peak->y - subImg->row0; + // + // Find the level nsigma below the peak that must separate the peak + // from any of its friends + // + assert (x >= 0 && x < subImg->numCols && y >= 0 && y < subImg->numRows); + const float stdev = sqrt(subWt->data.F32[y][x]); + float threshold = subImg->data.F32[y][x] - nsigma_delta*stdev; + if (threshold < min_threshold) { + threshold = min_threshold; + } + if (threshold > subImg->data.F32[y][x]) { + threshold = subImg->data.F32[y][x] - 10*FLT_EPSILON; + } + + psArray *peakFootprints = pmFindFootprints(subImg, threshold, 1); + psImage *idImg = pmSetFootprintArrayIDs(peakFootprints, true); + + const int peak_id = idImg->data.S32[y][x]; // the ID for the peak of interest + assert (peak_id != 0); + + int j; + for (j = 0; j < i; j++) { + const pmPeak *peak2 = fp->peaks->data[j]; + int x2 = peak2->x - subImg->col0; + int y2 = peak2->y - subImg->row0; + const int peak2_id = idImg->data.S32[y2][x2]; // the ID for some other peak + + if (peak2_id == peak_id) { // There's a brighter peak within the footprint above + ; // threshold; so cull our initial peak + (void)psArrayRemoveIndex(fp->peaks, i); + i--; // we moved everything down one + break; + } + } + if (j == i) { + j++; + } + + psFree(peakFootprints); + psFree(idImg); + } + + psFree((psImage *)subImg); + psFree((psImage *)subWt); + + return PS_ERR_NONE; +} + +/* + * Cull an entire psArray of pmFootprints + */ +psErrorCode +pmFootprintArrayCullPeaks(const psImage *img, // the image wherein lives the footprint + const psImage *weight, // corresponding variance image + psArray *footprints, // array of pmFootprints + const float nsigma_delta, // how many sigma above local background a peak + // needs to be to survive + const float min_threshold) { // minimum permitted coll height + for (int i = 0; i < footprints->n; i++) { + pmFootprint *fp = footprints->data[i]; + if (pmFootprintCullPeaks(img, weight, fp, nsigma_delta, min_threshold) != PS_ERR_NONE) { + return psError(PS_ERR_UNKNOWN, false, "Culling pmFootprint %d", fp->id); + } + } + + return PS_ERR_NONE; +} + +/************************************************************************************************************/ +/* + * Extract the peaks in a psArray of pmFootprints, returning a psArray of pmPeaks + */ +psArray *pmFootprintArrayToPeaks(const psArray *footprints) { + assert(footprints != NULL); + assert(footprints->n == 0 || pmIsFootprint(footprints->data[0])); + + int npeak = 0; + for (int i = 0; i < footprints->n; i++) { + const pmFootprint *fp = footprints->data[i]; + npeak += fp->peaks->n; + } + + psArray *peaks = psArrayAllocEmpty(npeak); + + for (int i = 0; i < footprints->n; i++) { + const pmFootprint *fp = footprints->data[i]; + for (int j = 0; j < fp->peaks->n; j++) { + psArrayAdd(peaks, 1, fp->peaks->data[j]); + } + } + + return peaks; +} Index: /trunk/psphot/src/pmFootprint.h =================================================================== --- /trunk/psphot/src/pmFootprint.h (revision 13007) +++ /trunk/psphot/src/pmFootprint.h (revision 13007) @@ -0,0 +1,43 @@ +#if !defined(PM_FOOTPRINT_H) +#define PM_FOOTPRINT_H + +typedef struct { + int y; + int x0, x1; +} pmSpan; + +pmSpan *pmSpanAlloc(int y, int x1, int x2); +bool pmIsSpan(const psPtr ptr); + +typedef struct { + const int id; // unique ID + int npix; // number of pixels in this pmFootprint + int nspan; // number of pmSpans in this pmFootprint + psArray *spans; // the pmSpans + psRegion bbox; // the pmFootprint's bounding box + psArray *peaks; // the peaks lying in this footprint + psRegion region; // A region describing the psImage the footprints live in +} pmFootprint; + +pmFootprint *pmFootprintAlloc(int nspan, const psImage *img); +bool pmIsFootprint(const psPtr ptr); + +int pmFootprintSetNpix(pmFootprint *fp); +void pmFootprintSetBBox(pmFootprint *fp); +psArray *pmFindFootprints(const psImage *img, const float threshold, const int npixMin); + +psArray *pmMergeFootprintArrays(const psArray *footprints1, const psArray *footprints2); + +psImage *pmSetFootprintArrayIDs(const psArray *footprints, + const bool relativeIDs); + +psErrorCode pmPeaksAssignToFootprints(psArray *footprints, const psArray *peaks); + +psErrorCode pmFootprintArrayCullPeaks(const psImage *img, const psImage *weight, psArray *footprints, + const float nsigma, const float threshold_min); +psErrorCode pmFootprintCullPeaks(const psImage *img, const psImage *weight, pmFootprint *fp, + const float nsigma, const float threshold_min); + +psArray *pmFootprintArrayToPeaks(const psArray *footprints); + +#endif Index: /trunk/psphot/src/psphot.h =================================================================== --- /trunk/psphot/src/psphot.h (revision 13006) +++ /trunk/psphot/src/psphot.h (revision 13007) @@ -30,5 +30,8 @@ // psphotReadout functions bool psphotImageMedian (pmConfig *config, pmFPAview *view); -psArray *psphotFindPeaks (pmReadout *readout, psMetadata *recipe, int pass); +psArray *psphotFindPeaks (pmReadout *readout, psMetadata *recipe, + const bool returnFootprints, const int pass); +#include "pmFootprint.h" +psErrorCode psphotCullPeaks(const pmReadout *readout, const psMetadata *recipe, psArray *footprints); psArray *psphotSourceStats (pmReadout *readout, psMetadata *recipe, psArray *allpeaks); bool psphotRoughClass (psArray *sources, psMetadata *recipe); @@ -48,6 +51,8 @@ // basic support functions void psphotModelGroupInit (void); -int psphotSortBySN (const void **a, const void **b); -int psphotSortByY (const void **a, const void **b); +int pmPeakSortBySN (const void **a, const void **b); +int pmPeakSortByY (const void **a, const void **b); +int pmSourceSortBySN (const void **a, const void **b); +int pmSourceSortByY (const void **a, const void **b); bool psphotGrowthCurve (pmReadout *readout, pmPSF *psf, bool ignore); bool psphotMaskReadout (pmReadout *readout, psMetadata *recipe); @@ -68,5 +73,5 @@ int psphotSaveImage (psMetadata *header, psImage *image, char *filename); bool psphotDumpConfig (pmConfig *config); -pmReadout *psphotSelectBackground (pmConfig *config, pmFPAview *view); +pmReadout *psphotSelectBackground (pmConfig *config, pmFPAview *view, const bool stdev); // PSF / DBL / EXT evaluation functions