IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 18847


Ignore:
Timestamp:
Aug 1, 2008, 8:39:03 AM (18 years ago)
Author:
Paul Price
Message:

Need to check return value of fgets to avoid compiler warning on some
systems.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/psphot/src/psphotExtendedSourceFits.c

    r18555 r18847  
    2424    // perform full extended source non-linear fits?
    2525    if (!psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_FITS")) {
    26         psLogMsg ("psphot", PS_LOG_INFO, "skipping extended source measurements\n");
    27         return true;
     26        psLogMsg ("psphot", PS_LOG_INFO, "skipping extended source measurements\n");
     27        return true;
    2828    }
    2929
     
    3131    psMetadata *models = psMetadataLookupMetadata (&status, recipe, "EXTENDED_SOURCE_MODELS");
    3232    if (!status) {
    33         psWarning ("extended source model fits requested but model model is missing (EXTENDED_SOURCE_MODELS)\n");
    34         return true;
     33        psWarning ("extended source model fits requested but model model is missing (EXTENDED_SOURCE_MODELS)\n");
     34        return true;
    3535    }
    3636    if (models->list->n == 0) {
    37         psWarning ("extended source model fits requested but no models are specified\n");
    38         return true;
     37        psWarning ("extended source model fits requested but no models are specified\n");
     38        return true;
    3939    }
    4040
     
    4545
    4646      if (item->type != PS_DATA_METADATA) {
    47         psAbort ("Invalid type for EXTENDED_SOURCE_MODEL entry %s, not a metadata folder", item->name);
    48         // XXX we could cull the bad entries or build a validated model folder
     47        psAbort ("Invalid type for EXTENDED_SOURCE_MODEL entry %s, not a metadata folder", item->name);
     48        // XXX we could cull the bad entries or build a validated model folder
    4949      }
    5050
     
    5555      int modelType = pmModelClassGetType (modelName);
    5656      if (modelType < 0) {
    57         psAbort ("Unknown model class for EXTENDED_SOURCE_MODEL entry %s: %s", item->name, modelName);
    58       } 
     57        psAbort ("Unknown model class for EXTENDED_SOURCE_MODEL entry %s: %s", item->name, modelName);
     58      }
    5959      psMetadataAddS32 (model, PS_LIST_TAIL, "MODEL_TYPE", PS_META_REPLACE, "", modelType);
    6060
     
    6262      char *SNword = psMetadataLookupStr (&status, model, "SNLIM");
    6363      if (!status) {
    64         psAbort("SNLIM not defined for extended source model %s\n", item->name);
     64        psAbort("SNLIM not defined for extended source model %s\n", item->name);
    6565      }
    6666      float SNlim = atof (SNword);
    6767      psMetadataAddF32 (model, PS_LIST_TAIL, "SNLIM_VALUE", PS_META_REPLACE, "", SNlim);
    68        
     68
    6969      // check on the PSF-Convolution status
    7070      char *convolvedWord = psMetadataLookupStr (&status, model, "PSF_CONVOLVED");
    7171      if (!status || (strcasecmp (convolvedWord, "true") && strcasecmp (convolvedWord, "false"))) {
    72         psAbort ("PSF_CONVOLVED entry invalid or missing for EXTENDED_SOURCE_MODEL entry %s", item->name);
    73       } 
     72        psAbort ("PSF_CONVOLVED entry invalid or missing for EXTENDED_SOURCE_MODEL entry %s", item->name);
     73      }
    7474      bool convolved = !strcasecmp (convolvedWord, "true");
    7575      psMetadataAddBool (model, PS_LIST_TAIL, "PSF_CONVOLVED_VALUE", PS_META_REPLACE, "", convolved);
     
    106106
    107107        // if model is NULL, we don't have a starting guess
    108         // XXX use the parameters guessed from moments
     108        // XXX use the parameters guessed from moments
    109109        // if (source->modelEXT == NULL) continue;
    110110
     
    115115        Next ++;
    116116
    117         // save the modelFlux here in case we need to subtract it (for failure)
    118         psImage *modelFluxStart = psMemIncrRefCounter (source->modelFlux);
    119 
    120         if (savePics) {
    121           psphotSaveImage (NULL, readout->image, "image.xp.fits");
    122         }
    123 
    124         // array to store the pointers to the model flux images while the models are being fitted
    125         psArray *modelFluxes = psArrayAllocEmpty (4);
    126    
    127         // allocate the array to store the model fits
    128         if (source->modelFits == NULL) {
    129             source->modelFits = psArrayAllocEmpty (4);
    130         }
    131 
    132         // loop here over the models chosen for each source (exclude by S/N)
    133         psMetadataIterator *iter = psMetadataIteratorAlloc (models, PS_LIST_HEAD, NULL);
    134         psMetadataItem *item = NULL;
    135         while ((item = psMetadataGetAndIncrement (iter)) != NULL) {
    136 
    137           // XXX this should have been forced above
    138           assert (item->type == PS_DATA_METADATA);
    139           psMetadata *model = (psMetadata *) item->data.md;
    140 
    141           // check the SNlim and skip model if source is too faint
    142           float SNlim = psMetadataLookupF32 (&status, model, "SNLIM_VALUE");
    143           assert (status);
    144 
    145           // limit selection to some SN limit
    146           assert (source->peak); // how can a source not have a peak?
    147           if (source->peak->SN < SNlim) continue;
    148 
    149           // check on the model type
    150           pmModelType modelType = psMetadataLookupS32 (&status, model, "MODEL_TYPE");
    151           assert (status);
    152 
    153           // check on the PSF-Convolution status
    154           bool convolved = psMetadataLookupBool (&status, model, "PSF_CONVOLVED_VALUE");
    155           assert (status);
    156 
    157           // XXX psTraceSetLevel ("psLib.math.psMinimizeLMChi2", 6);
    158           // XXX psTraceSetLevel ("psphot.psphotModelWithPSF_LMM", 6);
    159 
    160           // fit the model as convolved or not
    161           pmModel *modelFit = NULL;
    162           if (convolved) {
    163               modelFit = psphotPSFConvModel (readout, source, modelType, maskVal, markVal, psfSize);
    164               if (!modelFit) {
    165                   psTrace ("psphot", 5, "failed to fit psf-conv model for object at %f, %f", source->moments->x, source->moments->y);
    166                   continue;
    167               }
    168               psTrace ("psphot", 4, "fit psf-conv model for %f, %f : %s chisq = %f\n", source->moments->x, source->moments->y, pmModelClassGetName (modelFit->type), modelFit->chisq);
    169               Nconvolve ++;
    170               if (!(modelFit->flags & (PM_MODEL_STATUS_BADARGS | PM_MODEL_STATUS_NONCONVERGE | PM_MODEL_STATUS_OFFIMAGE))) {
    171                   NconvolvePass ++;
    172               }
    173           } else {
    174               psFree (source->modelFlux);
    175               source->modelFlux = NULL;
    176               modelFit = psphotFitEXT (readout, source, modelType, maskVal, markVal);
    177               if (!modelFit) {
    178                   psTrace ("psphot", 5, "failed to fit plain model for object at %f, %f", source->moments->x, source->moments->y);
    179                   continue;
    180               }
    181               pmSourceCacheModel (source, maskVal);
    182               psTrace ("psphot", 4, "fit plain model for %f, %f : %s chisq = %f\n", source->moments->x, source->moments->y, pmModelClassGetName (modelFit->type), modelFit->chisq);
    183               Nplain ++;
    184               if (!(modelFit->flags & (PM_MODEL_STATUS_BADARGS | PM_MODEL_STATUS_NONCONVERGE | PM_MODEL_STATUS_OFFIMAGE))) {
    185                   NplainPass ++;
    186               }
    187           }
    188 
    189           // save each of the model flux images and store the best
    190           psArrayAdd (modelFluxes, 4, source->modelFlux);
    191 
    192           // test for fit quality / result
    193           psArrayAdd (source->modelFits, 4, modelFit);
    194 
    195           psFree (modelFit);
    196         }
    197         psFree (iter);
    198 
    199         // evaluate the relative quality of the models, choose one
    200         float minChisq = NAN;
    201         int minModel = -1;
    202         for (int i = 0; i < source->modelFits->n; i++) {
    203             pmModel *model = source->modelFits->data[i];
    204            
    205             if (!(model->flags & PM_MODEL_STATUS_FITTED)) continue;
    206            
    207             if (model->flags & (PM_MODEL_STATUS_BADARGS)) continue;
    208             if (model->flags & (PM_MODEL_STATUS_NONCONVERGE)) continue;
    209             if (model->flags & (PM_MODEL_STATUS_OFFIMAGE)) continue;
    210 
    211             if ((minModel < 0) || (model->chisq < minChisq)) {
    212                 minChisq = model->chisq;
    213                 minModel = i;
    214             }
    215         }           
    216 
    217         if (minModel == -1) {
    218           // re-subtract the object, leave local sky
    219           psTrace ("psphot", 5, "failed to fit extended source model to object at %f, %f", source->moments->x, source->moments->y);
    220 
    221           // replace original model, subtract it
    222           psFree (source->modelFlux);
    223           source->modelFlux = modelFluxStart;
    224 
    225           pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
    226           source->mode |= PM_SOURCE_MODE_SUBTRACTED;
    227 
    228           psFree (modelFluxes);
    229 
    230           if (savePics) {
    231             psphotSaveImage (NULL, readout->image, "image.xp.fits");
    232             char key[10];
    233             fprintf (stdout, "continue? ");
    234             fgets (key, 8, stdin);
    235             if (key[0] == 'n') {
    236               savePics = false;
    237             }
    238           }
    239           continue;
    240         }       
    241 
    242         // save the best extended model in modelEXT
    243         psFree (source->modelEXT);
    244         source->modelEXT = psMemIncrRefCounter (source->modelFits->data[minModel]);
    245 
    246         // save the modelFlux for the best model
    247         psFree (source->modelFlux);
    248         source->modelFlux = psMemIncrRefCounter (modelFluxes->data[minModel]);
     117        // save the modelFlux here in case we need to subtract it (for failure)
     118        psImage *modelFluxStart = psMemIncrRefCounter (source->modelFlux);
     119
     120        if (savePics) {
     121          psphotSaveImage (NULL, readout->image, "image.xp.fits");
     122        }
     123
     124        // array to store the pointers to the model flux images while the models are being fitted
     125        psArray *modelFluxes = psArrayAllocEmpty (4);
     126
     127        // allocate the array to store the model fits
     128        if (source->modelFits == NULL) {
     129            source->modelFits = psArrayAllocEmpty (4);
     130        }
     131
     132        // loop here over the models chosen for each source (exclude by S/N)
     133        psMetadataIterator *iter = psMetadataIteratorAlloc (models, PS_LIST_HEAD, NULL);
     134        psMetadataItem *item = NULL;
     135        while ((item = psMetadataGetAndIncrement (iter)) != NULL) {
     136
     137          // XXX this should have been forced above
     138          assert (item->type == PS_DATA_METADATA);
     139          psMetadata *model = (psMetadata *) item->data.md;
     140
     141          // check the SNlim and skip model if source is too faint
     142          float SNlim = psMetadataLookupF32 (&status, model, "SNLIM_VALUE");
     143          assert (status);
     144
     145          // limit selection to some SN limit
     146          assert (source->peak); // how can a source not have a peak?
     147          if (source->peak->SN < SNlim) continue;
     148
     149          // check on the model type
     150          pmModelType modelType = psMetadataLookupS32 (&status, model, "MODEL_TYPE");
     151          assert (status);
     152
     153          // check on the PSF-Convolution status
     154          bool convolved = psMetadataLookupBool (&status, model, "PSF_CONVOLVED_VALUE");
     155          assert (status);
     156
     157          // XXX psTraceSetLevel ("psLib.math.psMinimizeLMChi2", 6);
     158          // XXX psTraceSetLevel ("psphot.psphotModelWithPSF_LMM", 6);
     159
     160          // fit the model as convolved or not
     161          pmModel *modelFit = NULL;
     162          if (convolved) {
     163              modelFit = psphotPSFConvModel (readout, source, modelType, maskVal, markVal, psfSize);
     164              if (!modelFit) {
     165                  psTrace ("psphot", 5, "failed to fit psf-conv model for object at %f, %f", source->moments->x, source->moments->y);
     166                  continue;
     167              }
     168              psTrace ("psphot", 4, "fit psf-conv model for %f, %f : %s chisq = %f\n", source->moments->x, source->moments->y, pmModelClassGetName (modelFit->type), modelFit->chisq);
     169              Nconvolve ++;
     170              if (!(modelFit->flags & (PM_MODEL_STATUS_BADARGS | PM_MODEL_STATUS_NONCONVERGE | PM_MODEL_STATUS_OFFIMAGE))) {
     171                  NconvolvePass ++;
     172              }
     173          } else {
     174              psFree (source->modelFlux);
     175              source->modelFlux = NULL;
     176              modelFit = psphotFitEXT (readout, source, modelType, maskVal, markVal);
     177              if (!modelFit) {
     178                  psTrace ("psphot", 5, "failed to fit plain model for object at %f, %f", source->moments->x, source->moments->y);
     179                  continue;
     180              }
     181              pmSourceCacheModel (source, maskVal);
     182              psTrace ("psphot", 4, "fit plain model for %f, %f : %s chisq = %f\n", source->moments->x, source->moments->y, pmModelClassGetName (modelFit->type), modelFit->chisq);
     183              Nplain ++;
     184              if (!(modelFit->flags & (PM_MODEL_STATUS_BADARGS | PM_MODEL_STATUS_NONCONVERGE | PM_MODEL_STATUS_OFFIMAGE))) {
     185                  NplainPass ++;
     186              }
     187          }
     188
     189          // save each of the model flux images and store the best
     190          psArrayAdd (modelFluxes, 4, source->modelFlux);
     191
     192          // test for fit quality / result
     193          psArrayAdd (source->modelFits, 4, modelFit);
     194
     195          psFree (modelFit);
     196        }
     197        psFree (iter);
     198
     199        // evaluate the relative quality of the models, choose one
     200        float minChisq = NAN;
     201        int minModel = -1;
     202        for (int i = 0; i < source->modelFits->n; i++) {
     203            pmModel *model = source->modelFits->data[i];
     204
     205            if (!(model->flags & PM_MODEL_STATUS_FITTED)) continue;
     206
     207            if (model->flags & (PM_MODEL_STATUS_BADARGS)) continue;
     208            if (model->flags & (PM_MODEL_STATUS_NONCONVERGE)) continue;
     209            if (model->flags & (PM_MODEL_STATUS_OFFIMAGE)) continue;
     210
     211            if ((minModel < 0) || (model->chisq < minChisq)) {
     212                minChisq = model->chisq;
     213                minModel = i;
     214            }
     215        }
     216
     217        if (minModel == -1) {
     218          // re-subtract the object, leave local sky
     219          psTrace ("psphot", 5, "failed to fit extended source model to object at %f, %f", source->moments->x, source->moments->y);
     220
     221          // replace original model, subtract it
     222          psFree (source->modelFlux);
     223          source->modelFlux = modelFluxStart;
     224
     225          pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
     226          source->mode |= PM_SOURCE_MODE_SUBTRACTED;
     227
     228          psFree (modelFluxes);
     229
     230          if (savePics) {
     231              psphotSaveImage (NULL, readout->image, "image.xp.fits");
     232              char key[10];
     233              fprintf (stdout, "continue? ");
     234              if (!fgets (key, 8, stdin)) {
     235                  psWarning("Couldn't read anything.");
     236              } else if (key[0] == 'n') {
     237                  savePics = false;
     238              }
     239          }
     240          continue;
     241        }
     242
     243        // save the best extended model in modelEXT
     244        psFree (source->modelEXT);
     245        source->modelEXT = psMemIncrRefCounter (source->modelFits->data[minModel]);
     246
     247        // save the modelFlux for the best model
     248        psFree (source->modelFlux);
     249        source->modelFlux = psMemIncrRefCounter (modelFluxes->data[minModel]);
    249250
    250251        // subtract the best fit from the object, leave local sky
     
    252253        source->mode |= PM_SOURCE_MODE_SUBTRACTED;
    253254
    254         // the initial model flux is no longer needed
    255         psFree (modelFluxStart);
    256         psFree (modelFluxes);
    257 
    258         psTrace ("psphot", 4, "best ext model for %f, %f : %s chisq = %f\n", source->moments->x, source->moments->y, pmModelClassGetName (source->modelEXT->type), source->modelEXT->chisq);
    259         psTrace ("psphot", 5, "extended source model for source at %7.1f, %7.1f", source->moments->x, source->moments->y);
    260 
    261         if (savePics) {
    262           psphotSaveImage (NULL, readout->image, "image.xm.fits");
    263           char key[10];
    264           fprintf (stdout, "continue? ");
    265           fgets (key, 8, stdin);
    266           if (key[0] == 'n') {
    267             savePics = false;
    268           }
    269         }
     255        // the initial model flux is no longer needed
     256        psFree (modelFluxStart);
     257        psFree (modelFluxes);
     258
     259        psTrace ("psphot", 4, "best ext model for %f, %f : %s chisq = %f\n", source->moments->x, source->moments->y, pmModelClassGetName (source->modelEXT->type), source->modelEXT->chisq);
     260        psTrace ("psphot", 5, "extended source model for source at %7.1f, %7.1f", source->moments->x, source->moments->y);
     261
     262        if (savePics) {
     263          psphotSaveImage (NULL, readout->image, "image.xm.fits");
     264          char key[10];
     265          fprintf (stdout, "continue? ");
     266          fgets (key, 8, stdin);
     267          if (key[0] == 'n') {
     268            savePics = false;
     269          }
     270        }
    270271    }
    271272
Note: See TracChangeset for help on using the changeset viewer.