IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 17702


Ignore:
Timestamp:
May 15, 2008, 12:01:36 PM (18 years ago)
Author:
eugene
Message:

working on getting the input image options right

Location:
branches/eam_branch_20080511/ppSim
Files:
11 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branch_20080511/ppSim/notes.txt

    r14657 r17702  
     1
     22008.05.14
     3
     4  ppSim needs some work to be a bit more realistic.  It currently does
     5  not do the correct thing for generating flats or shutter
     6  corrections.  It is only applying them to the background sky model.
     7  It is not applying them to the stellar flux.  Also, for adding
     8  sources on top of existing images, it needs to apply the Poisson
     9  stats to the generated stars and then add the pixel values to the
     10  input image.  currently, it is placing them directly in the input
     11  image.
     12
     13
    114
    215ppSim development work still needed:

  • branches/eam_branch_20080511/ppSim/src/Makefile.am

    r17672 r17702  
    2626        ppSimLoadSpots.c        \
    2727        ppSimPhotom.c           \
     28        ppSimMergeReadouts.c    \
    2829        ppSimMosaicChip.c       \
    2930        ppSimBadPixels.c

  • branches/eam_branch_20080511/ppSim/src/ppSim.h

    r17672 r17702  
    9595
    9696psVector *ppSimMakeBiassec (pmCell *cell, pmConfig *config);
    97 psVector *ppSimMakeBias (pmReadout *readout, pmConfig *config, const psRandom *rng) ;
     97psVector *ppSimMakeBias (bool *status, pmReadout *readout, pmConfig *config, const psRandom *rng) ;
    9898bool ppSimMakeDark (pmReadout *readout, pmConfig *config);
    9999bool ppSimMakeSky (pmReadout *readout, psImage *expCorr, ppSimType type, pmConfig *config);
    … …  
    174174float ppSimGetZeroPoint (psMetadata *recipe, char *filter);
    175175
     176bool ppSimMergeReadouts (pmReadout *outReadout, pmReadout *inReadout);
     177
    176178#endif

  • branches/eam_branch_20080511/ppSim/src/ppSimAddNoise.c

    r14657 r17702  
    4040    return true;
    4141}
    42 

  • branches/eam_branch_20080511/ppSim/src/ppSimAddOverscan.c

    r14657 r17702  
    44
    55    bool mdok;
     6
     7    // skip this step if we did not generate a bias level
     8    if (biasCols == NULL) return true;
    69
    710    psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSIM_RECIPE); // Recipe

  • branches/eam_branch_20080511/ppSim/src/ppSimArguments.c

    r17672 r17702  
    147147
    148148    char *typeStr = ppSimArgToRecipeStr (&status, options, "IMAGE.TYPE", arguments, "-type"); // Requested exposure type
    149     if (!status) {
     149    if (typeStr == NULL) {
    150150        psError(PS_ERR_IO, false, "An exposure type must be specified using -type");
    151151        psFree(arguments);
    … …  
    178178    ppSimArgToRecipeS32(&status, options, "BINNING",       arguments, "-bin");
    179179
     180    // XX don't do this: recipe is not yet fully realized (camera options not interpolated)
    180181    psMetadata *recipe = psMetadataLookupMetadata(&mdok, config->recipes, PPSIM_RECIPE); // Recipe
    181182

  • branches/eam_branch_20080511/ppSim/src/ppSimCreate.c

    r17672 r17702  
    1717            return NULL;
    1818        }
    19 
    20         simImage = true;
    21         fpa = pmFPAConstruct(config->camera); // FPA to contain the observation
    22         if (!fpa) {
    23             psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to construct an FPA from camera configuration.");
    24             return NULL;
    25         }
    2619    } else {
     20      // If an image is supplied, we still generate a fake image and merge them together downstream
     21      // (otherwise, we get the variance wrong).
    2722        simImage = false;
    2823        if (input->type != PM_FPA_FILE_IMAGE) {
    … …  
    3025            return NULL;
    3126        }
    32         fpa = psMemIncrRefCounter (input->fpa);
     27    }
     28
     29    // generate the fpa structure used by the output camera (determined from INPUT or specified)
     30    assert (config->camera);
     31    fpa = pmFPAConstruct(config->camera); // FPA to contain the observation
     32    if (!fpa) {
     33      psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to construct an FPA from camera configuration.");
     34      return NULL;
    3335    }
    3436
    … …  
    98100
    99101    // PPSIM.SOURCES carries the constructed, fake sources with their true parameters
    100     // XXX only invoke this code for OBJECT types of images
     102    // XXX only invoke this code for OBJECT types of images?
    101103    pmFPAfile *simSources = pmFPAfileDefineOutput (config, output->fpa, "PPSIM.SOURCES");
    102104    if (!simSources) {
    … …  
    106108    simSources->save = true;
    107109
    108     // if we have loaded an input image, we do not need to populate the fpa
    109     if (!simImage) {
     110    // if we have loaded an input image, we derive certain values from the image, if possible
     111    if (input) {
    110112        // we need to extract certain metadata from the image and populate the recipe.
    111113        // or else we need to set the fpa concepts based on the recipe options...
    … …  
    118120        float zp = ppSimGetZeroPoint (recipe, filter);
    119121        psMetadataAddF32(recipe, PS_LIST_TAIL, "ZEROPOINT", PS_META_REPLACE, "Photometric zeropoint", zp);
    120 
    121         return output;
    122122    }
    123123

  • branches/eam_branch_20080511/ppSim/src/ppSimLoop.c

    r17672 r17702  
    88bool ppSimLoop(pmConfig *config)
    99{
     10    bool status;
     11
    1012    PS_ASSERT_PTR_NON_NULL(config, PS_EXIT_PROG_ERROR);
    1113
    … …  
    1315    pmFPAfile *file = psMetadataLookupPtr(NULL, config->files, OUTPUT_FILE); // Output file
    1416    assert(file);
     17
     18    // if we have an input image, we need to add the synthetic data on top of it below
     19    // XXX we may potentially use this input file to skip missing elements
     20    pmFPAfile *input = psMetadataLookupPtr(NULL, config->files, "PPIMAGE.INPUT");
    1521
    1622    pmFPA *fpa = file->fpa;             // FPA for file
    … …  
    105111
    106112            psVector *biasCols = ppSimMakeBiassec (cell, config);
     113
     114            pmCell *inputCell = NULL;
     115            if (input) {
     116                inputCell = pmFPAviewThisCell (view, input->fpa);
     117            }
    107118
    108119            for (int i = 0; i < cell->readouts->n; i++) {
    … …  
    125136                }
    126137
    127                 psVector *biasRows = ppSimMakeBias (readout, config, rng);
     138                psVector *biasRows = ppSimMakeBias (&status, readout, config, rng);
     139                if (!status) ESCAPE (PS_ERR_UNKNOWN, "problem generating dark structure");
    128140                if (type == PPSIM_TYPE_BIAS) goto done;
    129141
    130                 ppSimMakeDark (readout, config);
     142                if (!ppSimMakeDark (readout, config)) ESCAPE (PS_ERR_UNKNOWN, "problem generating dark structure");
    131143                if (type == PPSIM_TYPE_DARK) goto done;
    132144
    133                 ppSimMakeSky (readout, expCorr, type, config);
     145                if (!ppSimMakeSky (readout, expCorr, type, config)) ESCAPE (PS_ERR_UNKNOWN, "problem generating sky background");
    134146                if (type == PPSIM_TYPE_FLAT) goto done;
    135147
    136148                if (type == PPSIM_TYPE_OBJECT) {
    137                     ppSimInsertStars (readout, expCorr, stars, config);
     149                    if (!ppSimInsertStars (readout, expCorr, stars, config)) ESCAPE (PS_ERR_UNKNOWN, "problem inserting stars");
    138150                }
    139151
    140152                if (type == PPSIM_TYPE_OBJECT) {
    141                     ppSimInsertGalaxies (readout, expCorr, galaxies, config);
     153                    if (!ppSimInsertGalaxies (readout, expCorr, galaxies, config)) ESCAPE (PS_ERR_UNKNOWN, "problem inserting galaxies");
    142154                }
    143155
    … …  
    145157
    146158            done:
    147                 ppSimAddNoise(readout->image, readout->weight, cell, config, rng);
    148                 ppSimSaturate(readout, config);
    149 
    150                 ppSimBadPixels(readout, config, rng);
    151 
    152                 ppSimAddOverscan (readout, config, biasCols, biasRows, rng);
     159                if (!ppSimAddNoise(readout->image, readout->weight, cell, config, rng)) ESCAPE (PS_ERR_UNKNOWN, "problem adding noise");
     160                if (!ppSimSaturate(readout, config)) ESCAPE (PS_ERR_UNKNOWN, "problem setting saturation levels");
     161
     162                if (!ppSimBadPixels(readout, config, rng)) ESCAPE (PS_ERR_UNKNOWN, "problem adding bad pixels");
     163
     164                if (!ppSimAddOverscan (readout, config, biasCols, biasRows, rng)) ESCAPE (PS_ERR_UNKNOWN, "problem adding overscan region");
    153165                psFree(biasRows);
    154166
    … …  
    156168                readout->parent->data_exists = true;
    157169                readout->parent->parent->data_exists = true;
     170
     171                if (inputCell && inputCell->readouts->data[i]) {
     172                    if (!ppSimMergeReadouts (readout, inputCell->readouts->data[i])) ESCAPE (PS_ERR_UNKNOWN, "problem merging input image with simulated image");
     173                }
    158174            }
    159175            psFree(biasCols);
    160176
    161             ppSimUpdateConceptsCell (cell, config);
     177            if (!ppSimUpdateConceptsCell (cell, config)) ESCAPE (PS_ERR_UNKNOWN, "problem updating cell concepts");
    162178
    163179            if (cell->hdu) {
    164                 ppSimInitHeader(config, NULL, NULL, cell);
     180                // XXX only do this if there is no INPUT image?
     181                if (!ppSimInitHeader(config, NULL, NULL, cell)) ESCAPE (PS_ERR_UNKNOWN, "problem setting output header");
    165182            }
    166183
    … …  
    175192
    176193        if (chip->hdu) {
    177             ppSimInitHeader(config, NULL, chip, NULL);
     194            // XXX only do this if there is no INPUT image
     195            if (!ppSimInitHeader(config, NULL, chip, NULL)) ESCAPE (PS_ERR_UNKNOWN, "problem setting output header");
    178196        }
    179197
    180198        // we perform photometry on the readouts of this chip in the output
    181         ppSimPhotom (config, view);
     199        if (!ppSimPhotom (config, view)) ESCAPE (PS_ERR_UNKNOWN, "problem performing photometry");
    182200
    183201        if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) {
    … …  
    185203            psFree(rng);
    186204            psFree(view);
    187             // return PS_EXIT_SYS_ERROR;
    188205            return false;
    189206        }
    … …  
    195212
    196213    if (fpa->hdu) {
    197         ppSimInitHeader(config, fpa, NULL, NULL);
     214        // XXX only do this if there is no INPUT image
     215        if (!ppSimInitHeader(config, fpa, NULL, NULL)) ESCAPE (PS_ERR_UNKNOWN, "problem setting output header");
    198216    }
    199217
    … …  
    202220        psFree(rng);
    203221        psFree(view);
    204         // return PS_EXIT_SYS_ERROR;
    205222        return false;
    206223    }

  • branches/eam_branch_20080511/ppSim/src/ppSimMakeBias.c

    r17672 r17702  
    11# include "ppSim.h"
    22
    3 psVector *ppSimMakeBias (pmReadout *readout, pmConfig *config, const psRandom *rng) {
     3psVector *ppSimMakeBias (bool *status, pmReadout *readout, pmConfig *config, const psRandom *rng) {
    44
    55    bool mdok;
     6
     7    if (status) *status = true;
    68
    79    pmCell *cell = readout->parent;
    … …  
    911    psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSIM_RECIPE); // Recipe
    1012
    11     bool bias = psMetadataLookupBool(&status, recipe, "BIAS"); // Generate a Bias?
    12     if (!bias) return true;
     13    bool bias = psMetadataLookupBool(&mdok, recipe, "BIAS"); // Generate a Bias?
     14    if (!bias) return NULL;
    1315
    1416    float biasLevel = psMetadataLookupF32(NULL, recipe, "BIAS.LEVEL"); // Bias level
    … …  
    2224        if (!mdok) {
    2325            psError(PS_ERR_BAD_PARAMETER_VALUE, false, "Unable to find READNOISE in recipe.");
     26            *status = false;
    2427            return NULL;
    2528        }

  • branches/eam_branch_20080511/ppSim/src/ppSimMakeDark.c

    r17672 r17702  
    1111    psMetadata *recipe = psMetadataLookupMetadata(&mdok, config->recipes, PPSIM_RECIPE); // Recipe
    1212
    13     bool dark = psMetadataLookupBool(&status, recipe, "DARK"); // Generate a DARK?
     13    bool dark = psMetadataLookupBool(&mdok, recipe, "DARK"); // Generate a DARK?
    1414    if (!dark) return true;
    1515
    16     float darkRate = psMetadataLookupF32(NULL, recipe, "DARK.RATE"); // Dark rate
    17     float expTime  = psMetadataLookupF32(NULL, recipe, "EXPTIME"); // Exposure time
     16    float darkRate = psMetadataLookupF32(&mdok, recipe, "DARK.RATE"); // Dark rate
     17    if (isnan(darkRate)) {
     18      psError(PS_ERR_UNKNOWN, false, "missing DARK.RATE\n");
     19      return false;
     20    }
     21
     22    float expTime  = psMetadataLookupF32(&mdok, recipe, "EXPTIME"); // Exposure time
     23    if (isnan(expTime)) {
     24      psError(PS_ERR_UNKNOWN, false, "missing EXPTIME\n");
     25      return false;
     26    }
     27
     28    psTrace("ppSim", 6, "darkRate: %f, expTime: %f\n", darkRate, expTime);
    1829
    1930    for (int y = 0; y < signal->numRows; y++) {

  • branches/eam_branch_20080511/ppSim/src/ppSimMakeSky.c

    r17672 r17702  
    11# include "ppSim.h"
    22
    3 // XXX add bounds to the inputs?
     3// this function sets the skyRate to a value for the night sky (SKY.RATE or SKY.MAGS) or for a
     4// flat-field image (FLAT.RATE).  Include a shutter correction and a scattered light source
     5
    46bool ppSimMakeSky (pmReadout *readout, psImage *expCorr, ppSimType type, pmConfig *config) {
    57
    … …  
    1618
    1719    bool sky = psMetadataLookupBool(&status, recipe, "SKY"); // Generate a SKY flux?
    18     if (!sky) return true;
    1920 
    20     float expTime     = psMetadataLookupF32(&status, recipe, "EXPTIME"); // Exposure time
    21     float flatSigma   = psMetadataLookupF32(&status, recipe, "FLAT.SIGMA"); // Flat-field coefficient
    22     float flatRate    = psMetadataLookupF32(&status, recipe, "FLAT.RATE"); // Flat-field rate
    23     float shutterTime = psMetadataLookupF32(&status, recipe, "SHUTTER.TIME"); // Shutter time
    24     float skyRate     = psMetadataLookupF32(&status, recipe, "SKY.RATE"); // Sky rate
     21    float expTime      = psMetadataLookupF32(&status, recipe, "EXPTIME"); // Exposure time
     22
     23    float flatSigma    = psMetadataLookupF32(&status, recipe, "FLAT.SIGMA"); // Flat-field coefficient
     24    float flatRate     = psMetadataLookupF32(&status, recipe, "FLAT.RATE"); // Flat-field rate
     25    float shutterTime  = psMetadataLookupF32(&status, recipe, "SHUTTER.TIME"); // Shutter time
     26    float scatterScale = psMetadataLookupF32(&status, recipe, "SCATTER.SCALE"); // scattered light term
     27    float skyRate      = psMetadataLookupF32(&status, recipe, "SKY.RATE"); // Sky rate
    2528    if (isnan(skyRate)) {
    26         float zp      = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); assert (status);
    27         float scale   = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE"); assert (status);
    28         float skyMags = psMetadataLookupF32(&status, recipe, "SKY.MAGS");  assert (status);
     29        float zp       = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); assert (status);
     30        float scale    = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE"); assert (status);
     31        float skyMags  = psMetadataLookupF32(&status, recipe, "SKY.MAGS");  assert (status);
    2932        skyRate = scale * scale * ppSimMagToFlux (skyMags, zp);
     33    }
     34    if (type == PPSIM_TYPE_FLAT) {
     35      skyRate = flatRate;
    3036    }
    3137
    … …  
    4450    // Size of FPA
    4551    psRegion *bounds = ppSimFPABounds (fpa);
     52    int dXfpa = bounds->x1 - bounds->x0;
     53    int dYfpa = bounds->y1 - bounds->y0;
    4654
    4755    // Correct chip offsets so that boresight is in the middle of the FPA
    48     x0Chip -= 0.5 * (bounds->x1 - bounds->x0);
    49     y0Chip -= 0.5 * (bounds->y1 - bounds->y0);
     56    x0Chip -= 0.5 * dXfpa;
     57    y0Chip -= 0.5 * dYfpa;
    5058
    5159    for (int y = 0; y < signal->numRows; y++) {
    … …  
    6270
    6371            // Gaussian flat-field over the FPA
    64             float flatValue = expf(-0.5 / PS_SQR(flatSigma) * (PS_SQR(yFPA) + PS_SQR(xFPA))) /
    65                 flatSigma / sqrtf(2.0 * M_PI);
     72            float flatValue = expf(-0.5 / PS_SQR(flatSigma) * (PS_SQR(yFPA) + PS_SQR(xFPA))) / flatSigma / sqrtf(2.0 * M_PI);
    6673
    67             if (type == PPSIM_TYPE_FLAT) {
    68                 float flatFlux = flatRate * flatValue * realExpTime; // Flux from flat-field
    69                 signal->data.F32[y][x] += flatFlux;
    70                 variance->data.F32[y][x] += flatFlux;
    71                 continue;
    72             }
     74            float scatterRate = 0.0;
    7375
    74             expCorr->data.F32[y][x] = realExpTime / expTime;
     76            if (sky) {
     77              // add a scattered light term to the flat-field images (no
     78              if (type == PPSIM_TYPE_FLAT) {
     79                float xF = 2.0*(xFPA / dXfpa) - 1.0;
     80                scatterRate = scatterScale * PS_SQR(xF);
     81              }
    7582
    76             // Sky background
    77             float skyFlux = skyRate * flatValue * realExpTime; // Flux from sky
    78             signal->data.F32[y][x] += skyFlux;
    79             variance->data.F32[y][x] += skyFlux;
     83              // Sky background
     84              float skyFlux = (skyRate * flatValue + scatterRate) * realExpTime; // Flux from sky
     85              signal->data.F32[y][x] += skyFlux;
     86              variance->data.F32[y][x] += skyFlux;
     87            }
     88
     89            // used later to modify the star and galaxy photometry
     90            if (expCorr) {
     91              expCorr->data.F32[y][x] = realExpTime / expTime;
     92            }
    8093
    8194            // TO DO: Add fringes
Note: See TracChangeset for help on using the changeset viewer.