IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Dec 17, 2010, 9:57:34 AM (16 years ago)
Author:
eugene
Message:

various modifications to psphotStack: add needed metadata to header; enable multiple PSF-matched images; work to speed up radial aperture photometry; radialAper is now an array (one element per PSF size); better feedback on extended fit skips; clarify processing: main analysis is on the input images (raw or convolved) while radial apertures is on the psf-matched images

File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20101205/psphot/src/psphotStackReadout.c

    r30027 r30101  
    11# include "psphotInternal.h"
    22
     3// we have 3 possible real filesets:
    34# define STACK_RAW "PSPHOT.STACK.INPUT.RAW"
    45# define STACK_CNV "PSPHOT.STACK.INPUT.CNV"
    5 # define STACK_OUT "PSPHOT.STACK.OUTPUT.IMAGE"
     6# define STACK_OUT "PSPHOT.STACK.OUTPUT.IMAGE"  /* the psf-matched image */
     7
     8// we have 3 files on which we operate:
     9// DET (detection image)       : nominally RAW (optionally CNV?)
     10// SRC (source analysis image) : nominally CNV (optionally RAW)
     11// OUT (psf-matched images)    : always OUT
    612
    713bool psphotStackReadout (pmConfig *config, const pmFPAview *view) {
    … …  
    2127    PS_ASSERT_PTR_NON_NULL (breakPt, false);
    2228
    23     BST = RAW : CNV;
     29    // XXX or do I set OUT to be a pmFPAfile pointing at the input of interest?
     30    bool useRaw = psMetadataLookupBool (NULL, recipe, "PSPHOT.STACK.USE.RAW");
     31    char *STACK_SRC = useRaw ? STACK_RAW : STACK_CNV;
     32    char *STACK_DET = STACK_RAW; // XXX optionally allow this to be CNV?
    2433
    2534    // we have 3 relevant files: RAW, CNV, OUT
    2635
    2736    // set the photcode for each image
    28     if (!psphotAddPhotcode (config, view, STACK_BST)) {
     37    if (!psphotAddPhotcode (config, view, STACK_SRC)) {
    2938        psError (PSPHOT_ERR_CONFIG, false, "trouble defining the photcode");
    3039        return false;
    … …  
    3342    // Generate the mask and weight images
    3443    // XXX this should be done before we perform the convolutions
    35     if (!psphotSetMaskAndVariance (config, view, STACK_RAW)) {
    36         return psphotReadoutCleanup (config, view, STACK_OUT);
     44    if (!psphotSetMaskAndVariance (config, view, STACK_DET)) {
     45        return psphotReadoutCleanup (config, view, STACK_SRC);
    3746    }
    3847    if (!strcasecmp (breakPt, "NOTHING")) {
    39         return psphotReadoutCleanup (config, view, STACK_OUT);
     48        return psphotReadoutCleanup (config, view, STACK_SRC);
    4049    }
    4150
    … …  
    4352    // XXX I think this is not defined correctly for an array of images.
    4453    // XXX probably need to subtract the model (same model?) for both RAW and OUT
    45     if (!psphotModelBackground (config, view, STACK_RAW)) {
    46         return psphotReadoutCleanup (config, view, STACK_OUT);
    47     }
    48     if (!psphotSubtractBackground (config, view, STACK_RAW)) {
    49         return psphotReadoutCleanup (config, view, STACK_OUT);
     54    if (!psphotModelBackground (config, view, STACK_DET)) {
     55        return psphotReadoutCleanup (config, view, STACK_SRC);
     56    }
     57    if (!psphotSubtractBackground (config, view, STACK_DET)) {
     58        return psphotReadoutCleanup (config, view, STACK_SRC);
    5059    }
    5160    if (!strcasecmp (breakPt, "BACKMDL")) {
    52         return psphotReadoutCleanup (config, view, STACK_OUT);
    53     }
    54 
    55     // load the psf model, if suppled.  FWHM_X,FWHM_Y,etc are determined and saved on
    56     // readout->analysis XXX this function currently only works with a single PSPHOT.INPUT
    57     if (!psphotLoadPSF (config, view, STACK_RAW)) {
    58         psError (PSPHOT_ERR_UNKNOWN, false, "error loading psf model");
    59         return psphotReadoutCleanup (config, view, STACK_OUT);
    60     }
    61 
    62     if (!psphotStackChisqImage(config, view, STACK_RAW, STACK_OUT)) {
     61        return psphotReadoutCleanup (config, view, STACK_SRC);
     62    }
     63
     64    if (!psphotStackChisqImage(config, view, STACK_DET, STACK_SRC)) {
    6365        psError (PSPHOT_ERR_UNKNOWN, false, "failure to generate chisq image");
    64         return psphotReadoutCleanup (config, view, STACK_OUT);
     66        return psphotReadoutCleanup (config, view, STACK_SRC);
    6567    }
    6668    if (!strcasecmp (breakPt, "CHISQ")) {
    67         return psphotReadoutCleanup (config, view, STACK_OUT);
     69        return psphotReadoutCleanup (config, view, STACK_SRC);
    6870    }
    6971
    7072    // find the detections (by peak and/or footprint) in the image.
    7173    // This finds the detections on Chisq image as well as the individuals
    72     if (!psphotFindDetections (config, view, STACK_RAW, true)) { // pass 1
     74    if (!psphotFindDetections (config, view, STACK_DET, true)) { // pass 1
    7375        // this only happens if we had an error in psphotFindDetections
    7476        psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis");
    75         return psphotReadoutCleanup (config, view, STACK_OUT);
    76     }
    77 
    78     // copy the detections from RAW to OUT
    79     if (!psphotCopySources (config, view, STACK_OUT, STACK_RAW)) {
    80         psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis");
    81         return psphotReadoutCleanup (config, view, STACK_OUT);
     77        return psphotReadoutCleanup (config, view, STACK_SRC);
     78    }
     79
     80    // copy the detections from DET to SRC
     81    if (strcmp(STACK_SRC, STACK_DET)) {
     82        if (!psphotCopySources (config, view, STACK_SRC, STACK_DET)) {
     83            psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis");
     84            return psphotReadoutCleanup (config, view, STACK_SRC);
     85        }
    8286    }
    8387
    8488    // construct sources and measure basic stats (saved on detections->newSources)
    85     // only run this on detections from the input images, not chisq image
    86     if (!psphotSourceStats (config, view, STACK_OUT, true)) { // pass 1
     89    if (!psphotSourceStats (config, view, STACK_SRC, true)) { // pass 1
    8790        psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources");
    88         return psphotReadoutCleanup (config, view, STACK_OUT);
     91        return psphotReadoutCleanup (config, view, STACK_SRC);
    8992    }
    9093
    9194    // generate the objects (object unify the sources from the different images)
    92     psArray *objects = psphotMatchSources (config, view, STACK_OUT);
     95    // XXX this could just match the detections for the chisq image, and not bother measuring the
     96    // source stats in that case...
     97    psArray *objects = psphotMatchSources (config, view, STACK_SRC);
    9398
    9499    // construct sources for the newly-generated sources (from other images)
    95     if (!psphotSourceStats (config, view, STACK_OUT, false)) { // pass 1
     100    if (!psphotSourceStats (config, view, STACK_SRC, false)) { // pass 1
    96101        psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources");
    97         return psphotReadoutCleanup (config, view, STACK_OUT);
     102        return psphotReadoutCleanup (config, view, STACK_SRC);
    98103    }
    99104
    … …  
    101106    // if (!psphotDeblendSatstars (config, view)) {
    102107    //     psError (PSPHOT_ERR_UNKNOWN, false, "failed on satstar deblend analysis");
    103     //     return psphotReadoutCleanup (config, view, STACK_OUT);
     108    //     return psphotReadoutCleanup (config, view, STACK_SRC);
    104109    // }
    105110
    … …  
    107112    // if (!psphotBasicDeblend (config, view)) {
    108113    //     psError (PSPHOT_ERR_UNKNOWN, false, "failed on deblend analysis");
    109     //     return psphotReadoutCleanup (config, view, STACK_OUT);
     114    //     return psphotReadoutCleanup (config, view, STACK_SRC);
    110115    // }
    111116
    112117    // classify sources based on moments, brightness
    113118    // only run this on detections from the input images, not chisq image
    114     if (!psphotRoughClass (config, view, STACK_OUT)) {
     119    if (!psphotRoughClass (config, view, STACK_SRC)) {
    115120        psError (PSPHOT_ERR_UNKNOWN, false, "failed to determine rough classifications");
    116         return psphotReadoutCleanup (config, view, STACK_OUT);
     121        return psphotReadoutCleanup (config, view, STACK_SRC);
    117122    }
    118123    // if we were not supplied a PSF model, determine the IQ stats here (detections->newSources)
    119124    // only run this on detections from the input images, not chisq image
    120     if (!psphotImageQuality (config, view, STACK_OUT)) { // pass 1
     125    if (!psphotImageQuality (config, view, STACK_SRC)) { // pass 1
    121126        psError (PSPHOT_ERR_UNKNOWN, false, "failed to measure image quality");
    122         return psphotReadoutCleanup (config, view, STACK_OUT);
     127        return psphotReadoutCleanup (config, view, STACK_SRC);
    123128    }
    124129    if (!strcasecmp (breakPt, "MOMENTS")) {
    125         return psphotReadoutCleanup (config, view, STACK_OUT);
     130        return psphotReadoutCleanup (config, view, STACK_SRC);
    126131    }
    127132
    128133    // use bright stellar objects to measure PSF if we were supplied a PSF for any input file,
    129134    // this step is skipped
    130     if (!psphotChoosePSF (config, view, STACK_OUT)) { // pass 1
     135    if (!psphotChoosePSF (config, view, STACK_SRC)) { // pass 1
    131136        psLogMsg ("psphot", 3, "failure to construct a psf model");
    132         return psphotReadoutCleanup (config, view, STACK_OUT);
     137        return psphotReadoutCleanup (config, view, STACK_SRC);
    133138    }
    134139    if (!strcasecmp (breakPt, "PSFMODEL")) {
    135         return psphotReadoutCleanup (config, view, STACK_OUT);
     140        return psphotReadoutCleanup (config, view, STACK_SRC);
    136141    }
    137142
    138143    // construct an initial model for each object, set the radius to fitRadius, set circular fit mask
    139     psphotGuessModels (config, view, STACK_OUT);
     144    psphotGuessModels (config, view, STACK_SRC);
    140145
    141146    // merge the newly selected sources into the existing list
    142147    // NOTE: merge OLD and NEW
    143     psphotMergeSources (config, view, STACK_OUT);
     148    psphotMergeSources (config, view, STACK_SRC);
    144149
    145150    // linear PSF fit to source peaks, subtract the models from the image (in PSF mask)
    … …  
    147152
    148153    // identify CRs and extended sources
    149     psphotSourceSize (config, view, STACK_OUT, TRUE);
     154    psphotSourceSize (config, view, STACK_SRC, TRUE);
    150155
    151156    // measure aperture photometry corrections
    152     if (!psphotApResid (config, view, STACK_OUT)) {
     157    if (!psphotApResid (config, view, STACK_SRC)) {
    153158        psFree (objects);
    154159        psLogMsg ("psphot", 3, "failed on psphotApResid");
    155         return psphotReadoutCleanup (config, view, STACK_OUT);
     160        return psphotReadoutCleanup (config, view, STACK_SRC);
    156161    }
    157162
    158163    psphotStackObjectsUnifyPosition (objects);
    159164
    160     // measure circular, radial apertures (objects sorted by S/N)
    161     psphotRadialAperturesByObject (config, objects, view, STACK_OUT);
    162 
    163165    // measure elliptical apertures, petrosians (objects sorted by S/N)
    164     psphotExtendedSourceAnalysisByObject (config, objects, view, STACK_OUT); // pass 1 (detections->allSources)
     166    psphotExtendedSourceAnalysisByObject (config, objects, view, STACK_SRC); // pass 1 (detections->allSources)
    165167
    166168    // measure non-linear extended source models (exponential, deVaucouleur, Sersic) (sources sorted by S/N)
    167     psphotExtendedSourceFits (config, view, STACK_OUT); // pass 1 (detections->allSources)
     169    psphotExtendedSourceFits (config, view, STACK_SRC); // pass 1 (detections->allSources)
    168170
    169171    // calculate source magnitudes
    170     psphotMagnitudes(config, view, STACK_OUT);
     172    psphotMagnitudes(config, view, STACK_SRC);
     173
     174    // copy the detections from SRC to OUT (for radial aperture photometry)
     175    if (!psphotCopySources (config, view, STACK_OUT, STACK_SRC)) {
     176        psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis");
     177        return psphotReadoutCleanup (config, view, STACK_SRC);
     178    }
     179    // XXX need to do something to reassign the source pixels here
     180
     181    bool smoothAgain = true;
     182    for (int nMatchedPSF = 0; smoothAgain; nMatchedPSF++) {
     183        // measure circular, radial apertures (objects sorted by S/N)
     184        psphotRadialAperturesByObject (config, objects, view, STACK_OUT, nMatchedPSF);
     185        psphotStackMatchPSFsNext(&smoothAgain, config, view, STACK_OUT, nMatchedPSF);
     186    }
    171187
    172188    if (0 && !psphotEfficiency(config, view, STACK_OUT)) {
    … …  
    179195
    180196    // replace background in residual image
    181     psphotSkyReplace (config, view, STACK_RAW);
     197    psphotSkyReplace (config, view, STACK_DET);
    182198
    183199    // drop the references to the image pixels held by each source
    184     psphotSourceFreePixels (config, view, STACK_OUT);
    185 
    186     // remove chisq image from config->file:PSPHOT.INPUT (why?)
    187     psphotStackRemoveChisqFromInputs(config, STACK_RAW);
     200    // psphotSourceFreePixels (config, view, STACK_OUT);
     201    psphotSourceFreePixels (config, view, STACK_SRC);
     202
     203    // remove chisq image from config->file:PSPHOT.INPUT
     204    psphotStackRemoveChisqFromInputs(config, STACK_DET);
     205    if (strcmp(STACK_SRC, STACK_DET)) {
     206        psphotStackRemoveChisqFromInputs(config, STACK_SRC);
     207    }
    188208
    189209    psFree (objects);
    190210
    191211    // create the exported-metadata and free local data
    192     return psphotReadoutCleanup (config, view, STACK_OUT);
     212    return psphotReadoutCleanup (config, view, STACK_SRC);
    193213}
    194214
    … …  
    258278 ******
    259279
    260  the above is all wrong:  first, we sohould be doing the full
     280 the above is all wrong:  first, we should be doing the full
    261281 morphology analysis (ExtendedAnalysis & ExtendedFits) on the CNV or
    262282 RAW image (as desired optionally), etc.
Note: See TracChangeset for help on using the changeset viewer.