IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Mar 4, 2011, 1:28:01 PM (15 years ago)
Author:
eugene
Message:

update the CELL.SATURATION value correctly based on inputs and recipe; plug some leaks

File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20110213/ppStack/src/ppStackConvolve.c

    r30620 r30801  
    44
    55// Update the value of a concept
    6 #define UPDATE_CONCEPT(SOURCE, NAME, VALUE) { \
    7     psMetadataItem *item = psMetadataLookup(SOURCE->concepts, NAME); \
    8     psAssert(item, "Concept should be present"); \
    9     psAssert(item->type == PS_TYPE_F32, "Concept should be F32"); \
    10     item->data.F32 = VALUE; \
    11 }
     6#define UPDATE_CONCEPT(SOURCE, NAME, VALUE) {                           \
     7        psMetadataItem *item = psMetadataLookup(SOURCE->concepts, NAME); \
     8        psAssert(item, "Concept should be present");                    \
     9        psAssert(item->type == PS_TYPE_F32, "Concept should be F32");   \
     10        item->data.F32 = VALUE;                                         \
     11    }
    1212
    1313bool ppStackConvolve(ppStackOptions *options, pmConfig *config)
    … …  
    4545    psVectorInit(renorms, NAN);
    4646
     47    psVector *satValues = psVectorAllocEmpty(num, PS_TYPE_F32); // Renormalisation values for variances
     48
    4749    psList *fpaList = psListAlloc(NULL); // List of input FPAs, for concept averaging
    4850    psList *cellList = psListAlloc(NULL); // List of input cells, for concept averaging
    … …  
    6870            psFree(cellList);
    6971            psFree(target);
     72            psFree(renorms);
     73            psFree(satValues);
    7074            return false;
    7175        }
    … …  
    8589            psFree(cellList);
    8690            psFree(target);
     91            psFree(renorms);
     92            psFree(satValues);
    8793            return false;
    8894        }
    … …  
    104110                psFree(cellList);
    105111                psFree(target);
     112                psFree(renorms);
     113                psFree(satValues);
    106114                return false;
    107115                // Non-fatal errors
    … …  
    118126                    psFree(cellList);
    119127                    psFree(target);
     128                    psFree(renorms);
     129                    psFree(satValues);
    120130                    return false;
    121131                }
    … …  
    164174            psFree(rng);
    165175            psFree(target);
     176            psFree(renorms);
     177            psFree(satValues);
    166178            return false;
    167179        }
    … …  
    175187            psFree(maskHeader);
    176188            psFree(target);
     189            psFree(renorms);
     190            psFree(satValues);
    177191            return false;
    178192        }
    … …  
    184198            psFree(rng);
    185199            psFree(target);
     200            psFree(renorms);
     201            psFree(satValues);
    186202            return false;
    187203        }
    … …  
    219235        if (options->matchZPs) {
    220236            // I think I need to take off the exposure time because we're going to set the new exposure time
     237            // Clarification: the zero point (ZP) in the header should be set such that:
     238            // M_app = m_inst + ZP + 2.5*log(exptime), where exptime in the output is sumExposure
    221239            psMetadataItem *zpItem = psMetadataLookup(inCell->parent->parent->concepts, "FPA.ZP");
    222240            zpItem->data.F32 += options->norm->data.F32[i] + 2.5*log10(options->sumExposure);
    … …  
    227245            expItem = psMetadataLookup(inCell->concepts, "CELL.EXPOSURE");
    228246            expItem->data.F32 = options->sumExposure;
     247
     248            // flux_out = flux_in * ten(-0.4*norm) -- save the individual saturation values
     249            psMetadataItem *satItem = psMetadataLookup(inCell->concepts, "CELL.EXPOSURE");
     250            satItem->data.F32 *= pow(10.0, -0.4*options->norm->data.F32[i]);
     251            psVectorAppend (satValues, satItem->data.F32);
    229252        }
    230253
    … …  
    235258            psFree(rng);
    236259            psFree(target);
     260            psFree(renorms);
     261            psFree(satValues);
    237262            return false;
    238263        }
    … …  
    255280        psFree(fpaList);
    256281        psFree(cellList);
     282        psFree(renorms);
     283        psFree(satValues);
    257284        return true;
    258285    }
    … …  
    279306        psFree(fpaList);
    280307        psFree(cellList);
    281 
     308       
     309        // The best guess for an output saturation value depends on the recipe.  If we have
     310        // 'safe' on, the we require at least 2 pixels to generate a valid output pixel.  In
     311        // this case, the best value for CELL.SATURATION is the 2nd highest value in the list.
     312        // If not, it should be the higest value in the list
     313        bool safe = psMetadataLookupBool(&mdok, recipe, "SAFE"); // Be safe when combining small numbers of pixels
     314        psVectorSortInPlace(satValues);
     315        float satBest = safe && satValues->n > 1 ? satValues->data.F32[1] : satValues->data.F32[0];
     316
     317        // UPDATE CELL.SATURATION here
    282318        UPDATE_CONCEPT(outFPA,  "FPA.EXPOSURE",  options->sumExposure);
    283319        UPDATE_CONCEPT(outCell, "CELL.EXPOSURE", options->sumExposure);
    284320        UPDATE_CONCEPT(outCell, "CELL.DARKTIME", NAN);
     321        UPDATE_CONCEPT(outCell, "CELL.SATURATION", satBest);
    285322        UPDATE_CONCEPT(outFPA,  "FPA.ZP",        options->zp);
    286323
    … …  
    288325        UPDATE_CONCEPT(unconvCell, "CELL.EXPOSURE", options->sumExposure);
    289326        UPDATE_CONCEPT(unconvCell, "CELL.DARKTIME", NAN);
     327        UPDATE_CONCEPT(unconvCell, "CELL.SATURATION", satBest);
    290328        UPDATE_CONCEPT(unconvFPA,  "FPA.ZP",        options->zp);
    291329
    … …  
    295333            double time = psTimeToMJD(fpaTime);
    296334            psMetadataAddF64(options->stats, PS_LIST_TAIL, "MJD_OBS", PS_META_DUPLICATE_OK,
    297                                  "Average MJD_OBS of inputs", time);
    298 
    299         }
    300     }
    301 
     335                             "Average MJD_OBS of inputs", time);
     336
     337        }
     338    }
     339
     340    // XXX EAM : this may be overly harsh -- or at least it would be if I (EAM) hadn't changed
     341    // the values of matchChi2 I modified pmSubtraction.c to fit a 2nd order polynomial to the
     342    // star chisq distribution (because of systematic errors in the model being matched).  This
     343    // fit forces the mean value to be 0.0.  Perhaps we can / should exclude images which have
     344    // an excessively high value for the rms?
    302345
    303346    // Reject images out-of-hand on the basis of their match chi^2
    … …  
    314357            psError(PPSTACK_ERR_PROG, false, "Unable to sort vector.");
    315358            psFree(values);
     359            psFree(renorms);
     360            psFree(satValues);
    316361            return false;
    317362        }
    … …  
    354399        psErrorClear();
    355400        psWarning("No good images survived convolution stage.");
     401        psFree(renorms);
     402        psFree(satValues);
    356403        return true;
    357404    }
    … …  
    364411    }
    365412    psFree(renorms);
     413    psFree(satValues);
    366414
    367415    return true;
Note: See TracChangeset for help on using the changeset viewer.