IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Jun 8, 2008, 4:03:31 PM (18 years ago)
Author:
eugene
Message:

update from changes on eam_branch_20080511 : adding photometry of fake sources, force photometry; major cleanups

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/ppSim/src/ppSimCreate.c

    r17915 r18011  
    11# include "ppSim.h"
    2 
    3 // XXX this function forces us to define the camera (on the command line) and format (via the
    4 // PPSIM.OUTPUT entry).  In this case, we need to set config->format,formatName based on these
    5 // values.  This will be a problem when we want to load an input image (in order to add fake
    6 // stars).  We will need to add some logic in ppSimArguments to distinguish the cases of 1)
    7 // input image and 2) specified camera
    82
    93pmFPAfile *ppSimCreate(pmConfig *config)
    104{
     5    PS_ASSERT_PTR_NON_NULL(config, NULL);
     6
    117    bool status;
    12     bool simImage = false;
    13     PS_ASSERT_PTR_NON_NULL(config, NULL);
    148    pmFPA *fpa = NULL;
    159
    1610    // the input image defines the camera.  if it is not supplied, the user must have
    1711    // supplied a camera and other metadata on the command line
    18     pmFPAfile *input = pmFPAfileDefineFromArgs (&status, config, "PPIMAGE.INPUT", "INPUT");
     12    pmFPAfile *input = pmFPAfileDefineFromArgs (&status, config, "PPSIM.INPUT", "INPUT");
    1913    if (!input) {
    20         // if we have not specified the camera already, we need to interpolate the recipes associated with this camera, and read other command-line recipes
    21         if (!pmConfigReadRecipes(config, PM_RECIPE_SOURCE_CL)) {
    22             psError(PS_ERR_IO, false, "Error merging recipes from camera config for %s", config->cameraName);
    23             return NULL;
    24         }
    25 
    26         simImage = true;
    27         fpa = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain the observation
    28         if (!fpa) {
    29             psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to construct an FPA from camera configuration.");
    30             return NULL;
    31         }
    32 
     14        // if we have not specified the camera already, we need to interpolate the recipes associated with this camera, and read other command-line recipes
     15        if (!pmConfigReadRecipes(config, PM_RECIPE_SOURCE_CL)) {
     16            psError(PS_ERR_IO, false, "Error merging recipes from camera config for %s", config->cameraName);
     17            return NULL;
     18        }
    3319    } else {
    34         simImage = false;
    35         if (input->type != PM_FPA_FILE_IMAGE) {
    36             psError(PS_ERR_IO, true, "PPIMAGE.INPUT is not of type IMAGE");
    37             return NULL;
    38         }
    39         fpa = input->fpa;
    40     }
    41 
    42     // define the output image file
    43     pmFPAfile *file = pmFPAfileDefineOutput(config, fpa, OUTPUT_FILE);
    44     if (!file) {
    45         psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to create output file from %s.  "
    46                 "Did you forget to specify the format?", OUTPUT_FILE);
    47         return NULL;
    48     }
    49     if (file->type != PM_FPA_FILE_IMAGE) {
    50         psError(PS_ERR_BAD_PARAMETER_TYPE, true, "%s type is not IMAGE", OUTPUT_FILE);
    51         psFree(fpa);
    52         psFree(file);
    53         return NULL;
    54     }
    55     file->save = true;
    56 
    57     config->format = psMemIncrRefCounter (file->format);
    58     config->formatName = psStringCopy (file->formatName);
    59 
    60     // have we supplied a psf model?
    61     if (psMetadataLookupPtr(NULL, config->arguments, "PSPHOT.PSF")) {
    62         bool status = false;
    63 
    64         // tie the psf file to the chipMosaic
    65         pmFPAfileBindFromArgs(&status, file, config, "PSPHOT.PSF.LOAD", "PSPHOT.PSF");
    66         if (!status) {
    67             psError(PS_ERR_UNKNOWN, false, "Failed to find/build PSPHOT.PSF.LOAD");
    68             psFree(fpa);
    69             psFree(file);
    70             return NULL;
    71         }
    72     }
    73 
    74     // XXX only invoke this code for OBJECT types of images
     20        // If an image is supplied, we still generate a fake image and merge them together downstream
     21        // (otherwise, we get the variance wrong).
     22        if (input->type != PM_FPA_FILE_IMAGE) {
     23            psError(PS_ERR_IO, true, "PPSIM.INPUT is not of type IMAGE");
     24            return NULL;
     25        }
     26    }
     27
     28    // generate the fpa structure used by the output camera (determined from INPUT or specified)
     29    assert (config->camera);
     30    fpa = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain the observation
     31    if (!fpa) {
     32        psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to construct an FPA from camera configuration.");
     33        return NULL;
     34    }
     35
     36    // define the output image file -- this is the basis for the ppSimLoop
     37    pmFPAfile *output = pmFPAfileDefineOutput(config, fpa, "PPSIM.OUTPUT");
     38    if (!output) {
     39        psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to create output file from PPSIM.OUTPUT. Did you forget to specify the format?");
     40        return NULL;
     41    }
     42    if (output->type != PM_FPA_FILE_IMAGE) {
     43        psError(PS_ERR_BAD_PARAMETER_TYPE, true, "PPSIM.OUTPUT type is not IMAGE");
     44        psFree(fpa);
     45        return NULL;
     46    }
     47    // XXX we should not require the output image to be written
     48    output->save = true;
     49
     50    config->format = psMemIncrRefCounter (output->format);
     51    config->formatName = psStringCopy (output->formatName);
     52
     53    // the recipe is now fully realized for the desired camera
     54    psMetadata *recipe = psMetadataLookupMetadata(&status, config->recipes, PPSIM_RECIPE); // Recipe
     55
     56    char *typeStr = psMetadataLookupStr(NULL, recipe, "IMAGE.TYPE"); // Type of image to simulate
     57    ppSimType type = ppSimTypeFromString (typeStr); // Type of image to simulate
     58
     59    if (type == PPSIM_TYPE_OBJECT) {
     60        // adjust the seeing by the scale
     61        float seeing = psMetadataLookupF32(&status, recipe, "SEEING");
     62        if (isnan(seeing)) {
     63            psError(PS_ERR_BAD_PARAMETER_TYPE, true, "seeing is not defined");
     64            psFree(fpa);
     65            return NULL;
     66        }
     67        float scale = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE");
     68        if (isnan(scale)) {
     69            psError(PS_ERR_BAD_PARAMETER_TYPE, true, "pixel scale is not defined");
     70            psFree(fpa);
     71            return NULL;
     72        }
     73        psMetadataAddF32(recipe, PS_LIST_TAIL, "SEEING", PS_META_REPLACE, "Seeing SIGMA (pixels)", seeing / 2.0 / sqrt(2.0 * log(2.0)) / scale);
     74
     75        // if we have been supplied an input image, but no ra & dec, use the input image values
     76        if (input) {
     77            float ra = psMetadataLookupF32(&status, recipe, "RA");
     78            if (isnan(ra)) {
     79                ra = psMetadataLookupF64(&status, input->fpa->concepts, "FPA.RA");
     80                psMetadataAddF32(recipe, PS_LIST_TAIL, "RA", PS_META_REPLACE, "ra (radians)", ra);
     81            }
     82            float dec = psMetadataLookupF32(&status, recipe, "DEC");
     83            if (isnan(dec)) {
     84                dec = psMetadataLookupF64(&status, input->fpa->concepts, "FPA.DEC");
     85                psMetadataAddF32(recipe, PS_LIST_TAIL, "DEC", PS_META_REPLACE, "dec (radians)", dec);
     86            }
     87        }
     88    }
     89
     90    if ((type == PPSIM_TYPE_OBJECT) || (type == PPSIM_TYPE_FLAT)) {
     91        // determine the zeropoint from the filter
     92        float zp = psMetadataLookupF32(&status, recipe, "ZEROPOINT");
     93        if (isnan(zp)) {
     94            char *filter = psMetadataLookupStr(&status, recipe, "FILTER");
     95            float zp = ppSimGetZeroPoint (recipe, filter);
     96            psMetadataAddF32(recipe, PS_LIST_TAIL, "ZEROPOINT", PS_META_REPLACE, "Photometric zeropoint", zp);
     97        }
     98    }
     99
     100    // For photometry, we operate on the chip-mosaicked image.  we create a copy of the mosaicked
     101    // image for psphot so we can write out a clean image
     102    bool doPhotom = psMetadataLookupBool(&status, recipe, "PHOTOM"); // Density of fakes
     103    if (doPhotom) {
     104
     105        // XXX at the moment, we can perform photometry on the fake sources and the forced
     106        // photometry positions, but not one or the other.  Also, we only support photometry in
     107        // the context of an image previously analysed by psphot.  Add support for:
     108        // * photometry of a pure fake image (requires peak detection and psf creation)
     109        // * photometry of forced source positions without fake image (this might work, it just
     110        // requires not generating any fake features).
     111
     112        if (!input) {
     113            psError(PS_ERR_UNKNOWN, false, "input image not found; currently required for photometry");
     114            return NULL;
     115        }
     116
     117        // we need a chip image if we perform photometry (is it necessary to build it if we don't use it?)
     118        pmFPAfile *fakeImage = pmFPAfileDefineChipMosaic(config, output->fpa, "PPSIM.FAKE.CHIP");
     119        if (!fakeImage) {
     120            psError(PS_ERR_IO, false, _("Unable to generate new file from PPSIM.FAKE.CHIP"));
     121            psFree(fpa);
     122            return NULL;
     123        }
     124        if (fakeImage->type != PM_FPA_FILE_IMAGE) {
     125            psError(PS_ERR_IO, true, "PPSIM.FAKE.CHIP is not of type IMAGE");
     126            psFree(fpa);
     127            return NULL;
     128        }
     129
     130        // we need a chip image if we perform photometry (is it necessary to build it if we don't use it?)
     131        pmFPAfile *forceImage = NULL;
     132        if (input) {
     133            forceImage = pmFPAfileDefineChipMosaic(config, input->fpa, "PPSIM.FORCE.CHIP");
     134            if (!forceImage) {
     135                psError(PS_ERR_IO, false, _("Unable to generate new file from PPSIM.FORCE.CHIP"));
     136                psFree(fpa);
     137                return NULL;
     138            }
     139            if (forceImage->type != PM_FPA_FILE_IMAGE) {
     140                psError(PS_ERR_IO, true, "PPSIM.FORCE.CHIP is not of type IMAGE");
     141                psFree(fpa);
     142                return NULL;
     143            }
     144        }
     145
     146        // define associated psphot input/output files
     147        if (!ppSimPhotomFiles (config, fakeImage, forceImage)) {
     148            psError(PSPHOT_ERR_CONFIG, false, "Trouble defining the additional input/output files for psphot");
     149            psFree(fpa);
     150            return NULL;
     151        }
     152    } else {
     153        // have we supplied a psf model?  this happens in ppSimPhotomFiles if we request a photometry
     154        // analysis.  however, even if we do not, a psf model may be used to generate the fake
     155        // sources.
     156        if (psMetadataLookupPtr(NULL, config->arguments, "PSPHOT.PSF")) {
     157            // tie the psf file to the chipMosaic
     158            pmFPAfileBindFromArgs(&status, output, config, "PSPHOT.PSF.LOAD", "PSPHOT.PSF");
     159            if (!status) {
     160                psError(PS_ERR_UNKNOWN, false, "Failed to find/build PSPHOT.PSF.LOAD");
     161                psFree(fpa);
     162                return NULL;
     163            }
     164        }
     165    }
     166
    75167    // PPSIM.SOURCES carries the constructed, fake sources with their true parameters
    76     pmFPAfile *simSources = pmFPAfileDefineOutput (config, file->fpa, "PPSIM.SOURCES");
     168    // XXX only invoke this code for OBJECT types of images?
     169    pmFPAfile *simSources = pmFPAfileDefineOutput (config, output->fpa, "PPSIM.SOURCES");
    77170    if (!simSources) {
    78         psError(PS_ERR_UNKNOWN, false, "Cannot find a rule for PPSIM.SOURCES");
    79         return false;
     171        psError(PS_ERR_UNKNOWN, false, "Cannot find a rule for PPSIM.SOURCES");
     172        return false;
    80173    }
    81174    simSources->save = true;
    82175
    83     // if we have loaded an input image, we do not need to populate the fpa
    84     if (!simImage) {
    85         return file;
    86     }
    87 
    88     pmFPALevel phuLevel = pmFPAPHULevel(file->format); // Level at which PHU goes
     176    // if we have loaded an input image, we derive certain values from the image, if possible
     177    if (input) {
     178        // we need to extract certain metadata from the image and populate the recipe.
     179        // or else we need to set the fpa concepts based on the recipe options...
     180
     181        psMetadata *recipe = psMetadataLookupMetadata(&status, config->recipes, PPSIM_RECIPE); // Recipe
     182
     183        ppSimArgToRecipeF32(&status, recipe, "EXPTIME", input->fpa->concepts, "FPA.EXPOSURE");
     184        char *filter = ppSimArgToRecipeStr(&status, recipe, "FILTER", input->fpa->concepts, "FPA.FILTER");
     185
     186        float zp = ppSimGetZeroPoint (recipe, filter);
     187        psMetadataAddF32(recipe, PS_LIST_TAIL, "ZEROPOINT", PS_META_REPLACE, "Photometric zeropoint", zp);
     188    }
     189
     190    pmFPALevel phuLevel = pmFPAPHULevel(output->format); // Level at which PHU goes
    89191
    90192    pmFPAview *view = pmFPAviewAlloc(0);// View for current level
    91193
    92194    if (phuLevel == PM_FPA_LEVEL_FPA) {
    93         if (!pmFPAAddSourceFromView(fpa, "Simulation", view, file->format)) {
    94             psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA.");
    95             psFree(fpa);
    96             psFree(view);
    97             psFree(file);
    98             return NULL;
    99         }
     195        if (!pmFPAAddSourceFromView(fpa, "Simulation", view, output->format)) {
     196            psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA.");
     197            psFree(fpa);
     198            psFree(view);
     199            return NULL;
     200        }
    100201    }
    101202
    102203    pmChip *chip;                       // Chip from FPA
    103204    while ((chip = pmFPAviewNextChip(view, fpa, 1))) {
    104         if (phuLevel == PM_FPA_LEVEL_CHIP) {
    105             if (!pmFPAAddSourceFromView(fpa, "Simulation", view, file->format)) {
    106                 psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA.");
    107                 psFree(fpa);
    108                 psFree(view);
    109                 psFree(file);
    110                 return NULL;
    111             }
    112         }
    113 
    114         pmCell *cell;                   // Cell from chip
    115         while ((cell = pmFPAviewNextCell(view, fpa, 1))) {
    116             if (phuLevel == PM_FPA_LEVEL_CELL) {
    117                 if (!pmFPAAddSourceFromView(fpa, "Simulation", view, file->format)) {
    118                     psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA.");
    119                     psFree(fpa);
    120                     psFree(view);
    121                     psFree(file);
    122                     return NULL;
    123                 }
    124             }
    125         }
     205        if (phuLevel == PM_FPA_LEVEL_CHIP) {
     206            if (!pmFPAAddSourceFromView(fpa, "Simulation", view, output->format)) {
     207                psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA.");
     208                psFree(fpa);
     209                psFree(view);
     210                return NULL;
     211            }
     212        }
     213
     214        pmCell *cell;                   // Cell from chip
     215        while ((cell = pmFPAviewNextCell(view, fpa, 1))) {
     216            if (phuLevel == PM_FPA_LEVEL_CELL) {
     217                if (!pmFPAAddSourceFromView(fpa, "Simulation", view, output->format)) {
     218                    psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA.");
     219                    psFree(fpa);
     220                    psFree(view);
     221                    return NULL;
     222                }
     223            }
     224        }
    126225    }
    127226
    … …  
    129228    psFree(view);
    130229
    131     return file;
     230    return output;
    132231}
Note: See TracChangeset for help on using the changeset viewer.