IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
May 18, 2010, 12:49:05 PM (16 years ago)
Author:
eugene
Message:

merging changes from trunk into branches/pap

Location:
branches/pap
Files:
20 edited
1 copied

Legend:

Unmodified
Added
Removed

  • branches/pap

  • branches/pap/psModules/src/astrom/pmAstrometryWCS.c

    r27461 r28003  
    972972    // FPA  -> TPA : identidy
    973973
     974    // NOTE: streaksremove's linearizeTransform function passes pointers to skeleton outFPA and outChip structs
     975    // Only outFPA->toSky is valid. No other memebers in these structs should be read. The resulting output transforms
     976    // are copied to inFPA and inChip by the caller.
     977
    974978    int nSamples = 10;  // 10 samples in each dimension
    975979

  • branches/pap/psModules/src/camera/pmFPAConstruct.c

    r22699 r28003  
    100100
    101101    return num;
    102 }
    103 
    104 // Get the name of a PHU chip or cell from the header
    105 static psString phuNameFromHeader(const char *name, // The name to lookup: "CELL.NAME" or "CHIP.NAME"
    106                                   const psMetadata *fileInfo, // FILE within the camera format description
    107                                   const psMetadata *header // Primary header
    108                                  )
    109 {
    110     assert(name && strlen(name) > 0);
    111     assert(fileInfo);
    112     assert(header);
    113 
    114     bool mdok = true;                   // Result of MD lookup
    115     psString keyword = psMetadataLookupStr(&mdok, fileInfo, name);
    116     if (!mdok || strlen(keyword) == 0) {
    117         return false;
    118     }
    119     psMetadataItem *resultItem = psMetadataLookup(header, keyword);
    120     if (!resultItem) {
    121         psError(PS_ERR_IO, true, "Unable to find %s in primary header to identify %s.\n", keyword, name);
    122         return NULL;
    123     }
    124     return psMetadataItemParseString(resultItem);
    125102}
    126103
    … …  
    11021079// It returns a view corresponding to the PHU
    11031080static pmFPAview *addSource(pmFPA *fpa,       // The FPA
    1104                             const char *fpaObs, // The desired FPA observation name
    11051081                            const pmFPAview *phuView, // The view corresponding to the PHU, or NULL
    11061082                            const psMetadata *header, // The PHU header, or NULL
    … …  
    11141090
    11151091    bool mdok;                          // Status of MD lookup
    1116 
    1117     // If FPA.OBS is already defined, new name must match it; otherwise, warn the user that something
    1118     // potentially bad is happening.
    1119     if (fpaObs && install) {
    1120         const char *currentName = psMetadataLookupStr(&mdok, fpa->concepts, "FPA.OBS"); // Current name
    1121         if (mdok && currentName && strlen(currentName) > 0 && strcmp(currentName, fpaObs) != 0) {
    1122             psWarning("FPA.OBS for new source (%s) doesn't match FPA.OBS for current fpa (%s).",
    1123                       fpaObs, currentName);
    1124         }
    1125         psMetadataAddStr(fpa->concepts, PS_LIST_HEAD, "FPA.OBS", PS_META_REPLACE, "Observation identifier",
    1126                          fpaObs);
    1127     } else if (!psMetadataLookup(fpa->concepts, "FPA.OBS")) {
    1128         // Make sure there is an FPA.OBS
    1129         psMetadataAddStr(fpa->concepts, PS_LIST_HEAD, "FPA.OBS", 0, "Observation identifier", "UNKNOWN");
    1130     }
    11311092
    11321093    psMetadata *fileInfo = psMetadataLookupMetadata(&mdok, format, "FILE"); // The file information
    … …  
    13501311}
    13511312
    1352 bool pmFPAAddSourceFromFormat(pmFPA *fpa, const char *fpaObs, const psMetadata *format)
     1313bool pmFPAAddSourceFromFormat(pmFPA *fpa, const psMetadata *format)
    13531314{
    13541315    PS_ASSERT_PTR_NON_NULL(fpa, false);
    … …  
    13591320    pmFPAview *view = pmFPAviewAlloc(0);// View for current level
    13601321    if (phuLevel == PM_FPA_LEVEL_FPA) {
    1361         if (!pmFPAAddSourceFromView(fpa, fpaObs, view, format)) {
     1322        if (!pmFPAAddSourceFromView(fpa, view, format)) {
    13621323            psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA.");
    13631324            psFree(view);
    … …  
    13681329        while ((chip = pmFPAviewNextChip(view, fpa, 1))) {
    13691330            if (phuLevel == PM_FPA_LEVEL_CHIP) {
    1370                 if (!pmFPAAddSourceFromView(fpa, fpaObs, view, format)) {
     1331                if (!pmFPAAddSourceFromView(fpa, view, format)) {
    13711332                    psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA.");
    13721333                    psFree(view);
    … …  
    13771338                while ((cell = pmFPAviewNextCell(view, fpa, 1))) {
    13781339                    if (phuLevel == PM_FPA_LEVEL_CELL) {
    1379                         if (!pmFPAAddSourceFromView(fpa, fpaObs, view, format)) {
     1340                        if (!pmFPAAddSourceFromView(fpa, view, format)) {
    13801341                            psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA.");
    13811342                            psFree(view);
    … …  
    13921353}
    13931354
    1394 bool pmFPAAddSourceFromView(pmFPA *fpa, const char *fpaObs, const pmFPAview *phuView,
     1355bool pmFPAAddSourceFromView(pmFPA *fpa, const pmFPAview *phuView,
    13951356                            const psMetadata *format)
    13961357{
    … …  
    13991360    PS_ASSERT_PTR_NON_NULL(format, false);
    14001361
    1401     pmFPAview *view = addSource(fpa, fpaObs, phuView, NULL, format, true);
     1362    pmFPAview *view = addSource(fpa, phuView, NULL, format, true);
    14021363    bool status = (view != NULL);
    14031364    psFree(view);
    … …  
    14181379    }
    14191380
    1420     // Check the name of the FPA
    1421     psString fpaObs = phuNameFromHeader("FPA.OBS", fileInfo, phu); // New observation name for the FPA
    1422     if (!fpaObs || strlen(fpaObs) == 0) {
    1423         psWarning("Unable to determine FPA.OBS: check for FPA.OBS in FILE in camera format");
    1424     }
    1425 
    1426     pmFPAview *view = addSource(fpa, fpaObs, NULL, phu, format, true); // View of PHU, to return
    1427     psFree(fpaObs);
     1381    pmFPAview *view = addSource(fpa, NULL, phu, format, true); // View of PHU, to return
    14281382
    14291383    return view;
    … …  
    14371391    PS_ASSERT_PTR_NON_NULL(format, NULL);
    14381392
    1439     return addSource(fpa, NULL, NULL, phu, format, false);
     1393    return addSource(fpa, NULL, phu, format, false);
    14401394}
    14411395

  • branches/pap/psModules/src/camera/pmFPAConstruct.h

    r17911 r28003  
    2929/// This is suitable for generating an output FPA given the desired format.
    3030bool pmFPAAddSourceFromFormat(pmFPA *fpa, ///< The FPA
    31                               const char *fpaObs, ///< FPA.NAME for the source
    3231                              const psMetadata *format ///< Format of file
    3332    );
    … …  
    3837/// configuration is required in order to describe how the FPA is laid out in terms of disk files.
    3938bool pmFPAAddSourceFromView(pmFPA *fpa,   ///< The FPA
    40                             const char *fpaObs, ///< FPA.NAME for the source
    4139                            const pmFPAview *phuView, ///< The view, corresponding to the PHU
    4240                            const psMetadata *format ///< Format of file

  • branches/pap/psModules/src/camera/pmFPAWrite.c

    r25881 r28003  
    165165    // generate the HDU, but only copies the structure.
    166166    if (!blank && !hdu->blankPHU && !*imageArray) {
    167         if (!pmHDUGenerateForCell(cell)) {
    168             psError(PS_ERR_UNKNOWN, false, "Unable to generate HDU for cell --- likely programming error.");
    169             return false;
    170         }
    171         if (!*imageArray) {
    172             if (type == FPA_WRITE_TYPE_IMAGE) {
    173                 psError(PS_ERR_UNKNOWN, false, "Expected to write an image, but it is missing...programming error?.");
    174                 return false;
    175             }
    176             if (type == FPA_WRITE_TYPE_MASK) {
    177                 psWarning("No mask image for this cell; skipping");
    178             }
    179             if (type == FPA_WRITE_TYPE_VARIANCE) {
    180                 psWarning("No variance image for this cell; skipping");
    181             }
    182             return true;
    183         }
     167        if (!pmHDUGenerateForCell(cell)) {
     168            psError(PS_ERR_UNKNOWN, false, "Unable to generate HDU for cell --- likely programming error.");
     169            return false;
     170        }
     171        if (!*imageArray) {
     172            if (type == FPA_WRITE_TYPE_IMAGE) {
     173                psError(PS_ERR_UNKNOWN, false, "Expected to write an image, but it is missing...programming error?.");
     174                return false;
     175            }
     176            if (type == FPA_WRITE_TYPE_MASK) {
     177                psWarning("No mask image for this cell; skipping");
     178            }
     179            if (type == FPA_WRITE_TYPE_VARIANCE) {
     180                psWarning("No variance image for this cell; skipping");
     181            }
     182            return true;
     183        }
    184184    }
    185185
    … …  
    307307
    308308        if (writeBlank || writeImage) {
    309             if (!pmConceptsWriteFPA(fpa, true, config)) {
     309            if (!pmConceptsWriteFPA(fpa, true, config)) {
    310310                psError(PS_ERR_IO, false, "Unable to write concepts for FPA.\n");
    311311                return false;
    … …  
    341341//////////////////////////////////////////////////////////////////////////////////////////////////////////////
    342342
    343 // Update the FPA.OBS, CHIP.NAME and CELL.NAME in the FITS header, if required
     343// Update the CHIP.NAME and CELL.NAME in the FITS header, if required
    344344bool pmFPAUpdateNames(pmFPA *fpa, pmChip *chip, pmCell *cell, psS64 imageId, psS64 sourceId)
    345345{
    … …  
    375375        psError(PS_ERR_UNKNOWN, true, "Unable to find FILE information in camera format.\n");
    376376        return false;
    377     }
    378     if (fpa) {
    379         const char *fpaObsHdr = psMetadataLookupStr(&mdok, fileData, "FPA.OBS");
    380         if (mdok && fpaObsHdr && strlen(fpaObsHdr) > 0) {
    381             const char *fpaObs = psMetadataLookupStr(NULL, fpa->concepts, "FPA.OBS");
    382             psMetadataAddStr(hduHigh->header, PS_LIST_TAIL, fpaObsHdr, PS_META_REPLACE,
    383                              "Observation identifier", fpaObs);
    384         }
    385377    }
    386378

  • branches/pap/psModules/src/camera/pmFPAfileFitsIO.c

    r25878 r28003  
    7575    }
    7676
    77     pmFPA *nameSource = file->src; // Source of FPA.OBS
    78     if (!nameSource) {
    79         nameSource = file->fpa;
    80     }
    81     bool mdok;                  // Status of MD lookup
    82     const char *fpaObs = psMetadataLookupStr(&mdok, nameSource->concepts, "FPA.OBS"); // Observation id
    83 
    8477    pmFPA *copy = pmFPAConstruct(file->camera, file->cameraName);  // FPA to return
    85     if (!pmFPAAddSourceFromView(copy, fpaObs, phuView, file->format)) {
     78    if (!pmFPAAddSourceFromView(copy, phuView, file->format)) {
    8679        psError(PS_ERR_UNKNOWN, false, "Unable to insert HDU into FPA for writing.\n");
    8780        psFree(copy);

  • branches/pap/psModules/src/camera/pmFPAfileIO.c

    r26893 r28003  
    286286            }
    287287
    288             pmFPA *nameSource = file->src; // Source of FPA.OBS
    289             if (!nameSource) {
    290                 nameSource = file->fpa;
    291             }
    292             bool mdok;                  // Status of MD lookup
    293             const char *fpaObs = psMetadataLookupStr(&mdok, nameSource->concepts, "FPA.OBS"); // Obs. id
    294 
    295             pmFPAAddSourceFromView(file->fpa, fpaObs, view, format);
     288            pmFPAAddSourceFromView(file->fpa, view, format);
    296289            psTrace ("psModules.camera", 5, "created fpa data elements for %s (%s) (%d:%d:%d)\n",
    297290                     file->name, file->name, view->chip, view->cell, view->readout);

  • branches/pap/psModules/src/detrend/pmShutterCorrection.c

    r23989 r28003  
    10151015
    10161016        pmShutterCorrection *corr = pmShutterCorrectionFullFit(newtimes, newcounts, newerrors, guess); // The actual correction
    1017         psTrace("psModules.detrend", 5, "Shutter correction fit: scale: %f, offset: %f, offref: %f\n", corr->scale, corr->offset, corr->offref);
    1018 
    1019         if (corr && isfinite(corr->offref) && corr->valid) {
    1020             psTrace("psModules.detrend", 5, "Sample reference value: %f\n", corr->offref);
    1021             meanRef += corr->offref;
    1022             numGood++;
    1023         }
     1017
     1018        if (corr) {
     1019            psTrace("psModules.detrend", 5, "Shutter correction fit: scale: %f, offset: %f, offref: %f\n", corr->scale, corr->offset, corr->offref);
     1020            if (isfinite(corr->offref) && corr->valid) {
     1021                psTrace("psModules.detrend", 5, "Sample reference value: %f\n", corr->offref);
     1022                meanRef += corr->offref;
     1023                numGood++;
     1024            }
     1025        } else {
     1026            psTrace("psModules.detrend", 5, "failed Shutter correction fit\n");
     1027        }
    10241028
    10251029        psFree(corr);

  • branches/pap/psModules/src/extras/Makefile.am

    r24880 r28003  
    1717        pmKapaPlots.h \
    1818        pmVisual.h \
     19        ippDiffMode.h \
    1920        ippStages.h
    2021

  • branches/pap/psModules/src/extras/ippStages.h

    r24880 r28003  
    1 /* @file psVisual.h
     1/* @file ippStages.h
    22 * @brief some macro defintions for the stages of the pipeline
    33 * @author Bill Sweeney, IfA

  • branches/pap/psModules/src/objects/pmPSFtryModel.c

    r27330 r28003  
    107107    // as we loop over orders, we need to refer to the initial selection, but we modify the
    108108    // option values to match the current guess: save the max values here:
    109     int Nx = options->psfTrendNx;
    110     int Ny = options->psfTrendNy;
     109    int Nx = (options->psfTrendMode == PM_TREND_MAP) ? options->psfTrendNx : options->psfTrendNx + 1;
     110    int Ny = (options->psfTrendMode == PM_TREND_MAP) ? options->psfTrendNy : options->psfTrendNy + 1;
    111111    for (int i = orderMin; i <= orderMax; i++) {
    112112

  • branches/pap/psModules/src/objects/pmSourceExtendedPars.c

    r25754 r28003  
    3636#include "pmSourceExtendedPars.h"
    3737
    38 // *** pmSourceRadialProfile describes the radial profile of a source in elliptical contours, and
    39 // intermediate data used to measure the profile
     38// pmSourceRadialFlux carries the raw radial flux information, including angular bins
     39static void pmSourceRadialFluxFree(pmSourceRadialFlux *flux)
     40{
     41    if (!flux) return;
     42    psFree(flux->radii);
     43    psFree(flux->fluxes);
     44    psFree(flux->theta);
     45    psFree(flux->isophotalRadii);
     46}
     47
     48pmSourceRadialFlux *pmSourceRadialFluxAlloc()
     49{
     50    pmSourceRadialFlux *flux = (pmSourceRadialFlux *)psAlloc(sizeof(pmSourceRadialFlux));
     51    psMemSetDeallocator(flux, (psFreeFunc) pmSourceRadialFluxFree);
     52
     53    flux->radii = NULL;
     54    flux->fluxes = NULL;
     55    flux->theta = NULL;
     56    flux->isophotalRadii = NULL;
     57
     58    return flux;
     59}
     60
     61bool psMemCheckSourceRadialFlux(psPtr ptr)
     62{
     63    PS_ASSERT_PTR(ptr, false);
     64    return ( psMemGetDeallocator(ptr) == (psFreeFunc) pmSourceRadialFluxFree);
     65}
     66
     67// pmSourceEllipticalFlux carries the elliptical renormalized radial flux info
     68static void pmSourceEllipticalFluxFree(pmSourceEllipticalFlux *flux)
     69{
     70    if (!flux) return;
     71    psFree(flux->radiusElliptical);
     72    psFree(flux->fluxElliptical);
     73}
     74
     75pmSourceEllipticalFlux *pmSourceEllipticalFluxAlloc()
     76{
     77    pmSourceEllipticalFlux *flux = (pmSourceEllipticalFlux *)psAlloc(sizeof(pmSourceEllipticalFlux));
     78    psMemSetDeallocator(flux, (psFreeFunc) pmSourceEllipticalFluxFree);
     79
     80    flux->radiusElliptical = NULL;
     81    flux->fluxElliptical = NULL;
     82
     83    return flux;
     84}
     85
     86bool psMemCheckSourceEllipticalFlux(psPtr ptr)
     87{
     88    PS_ASSERT_PTR(ptr, false);
     89    return ( psMemGetDeallocator(ptr) == (psFreeFunc) pmSourceEllipticalFluxFree);
     90}
     91
     92// pmSourceRadialProfile defines flux information in radial bins
    4093static void pmSourceRadialProfileFree(pmSourceRadialProfile *profile)
    4194{
    42     if (!profile) {
    43         return;
    44     }
    45     psFree(profile->radii);
    46     psFree(profile->fluxes);
    47     psFree(profile->theta);
    48     psFree(profile->isophotalRadii);
    49 
    50     psFree(profile->radiusElliptical);
    51     psFree(profile->fluxElliptical);
    52 
     95    if (!profile) return;
    5396    psFree(profile->binSB);
    5497    psFree(profile->binSBstdev);
    5598    psFree(profile->binSBerror);
    56 
     99    psFree(profile->binSum);
     100    psFree(profile->binFill);
    57101    psFree(profile->radialBins);
    58102    psFree(profile->area);
    … …  
    63107    pmSourceRadialProfile *profile = (pmSourceRadialProfile *)psAlloc(sizeof(pmSourceRadialProfile));
    64108    psMemSetDeallocator(profile, (psFreeFunc) pmSourceRadialProfileFree);
    65 
    66     profile->radii = NULL;
    67     profile->fluxes = NULL;
    68     profile->theta = NULL;
    69     profile->isophotalRadii = NULL;
    70 
    71     profile->radiusElliptical = NULL;
    72     profile->fluxElliptical = NULL;
    73109
    74110    profile->binSB = NULL;
    75111    profile->binSBstdev = NULL;
    76112    profile->binSBerror = NULL;
    77 
     113    profile->binSum = NULL;
     114    profile->binFill = NULL;
    78115    profile->radialBins = NULL;
    79116    profile->area = NULL;
    80 
    81117    return profile;
    82118}
    … …  
    88124}
    89125
    90 
    91 // *** pmSourceRadialProfileFreeVectors frees the intermediate data values
     126# if (0)
     127// pmSourceRadialProfileFreeVectors frees the intermediate data values
    92128bool pmSourceRadialProfileFreeVectors(pmSourceRadialProfile *profile) {
    93129
    … …  
    124160    return true;
    125161}
     162# endif
    126163
    127164// *** pmSourceRadialProfileSortPair is a utility function for sorting a pair of vectors
    … …  
    150187    if (!pars) return;
    151188
    152     psFree(pars->profile);
    153     psFree(pars->petrosian_50);
    154     psFree(pars->petrosian_80);
     189    psFree(pars->radFlux);
     190    psFree(pars->ellipticalFlux);
     191    psFree(pars->radProfile);
     192    psFree(pars->petProfile);
    155193    return;
    156194}
    … …  
    160198    psMemSetDeallocator(pars, (psFreeFunc) pmSourceExtendedParsFree);
    161199
    162     pars->profile = NULL;
    163     pars->petrosian_50 = NULL;
    164     pars->petrosian_80 = NULL;
    165 
     200    pars->radFlux = NULL;
     201    pars->ellipticalFlux = NULL;
     202    pars->radProfile = NULL;
     203    pars->petProfile = NULL;
     204
     205    pars->petrosianFlux = NAN;
     206    pars->petrosianFluxErr = NAN;
     207    pars->petrosianRadius = NAN;
     208    pars->petrosianRadiusErr = NAN;
     209    pars->petrosianR90 = NAN;
     210    pars->petrosianR90Err = NAN;
     211    pars->petrosianR50 = NAN;
     212    pars->petrosianR50Err = NAN;
    166213    return pars;
    167214}

  • branches/pap/psModules/src/objects/pmSourceExtendedPars.h

    r25754 r28003  
    1919    psVector *theta;                    // angles corresponding to above radial profiles
    2020    psVector *isophotalRadii;           // isophotal radius for the above angles
     21} pmSourceRadialFlux;
    2122
     23typedef struct {
    2224    psVector *radiusElliptical;         // normalized radial coordinates for all relevant pixels
    2325    psVector *fluxElliptical;           // flux for the above radial coordinates
     26} pmSourceEllipticalFlux;
    2427
     28typedef struct {
    2529    psVector *binSB;                    // mean surface brightness within radial bins
    2630    psVector *binSBstdev;               // scatter of mean surface brightness within radial bins
    2731    psVector *binSBerror;               // formal error on mean surface brightness within radial bins
    28 
    29     psVector *radialBins;               // radii corresponding to above binnedBlux
     32    psVector *binSum;                   // sum of flux within radial bins
     33    psVector *binFill;                  // fraction of area actually lit
     34    psVector *radialBins;               // radii corresponding to above binnedFlux
    3035    psVector *area;                     // differential area of the non-overlapping radial bins
    31 
    32     psEllipseAxes axes;                 // shape of elliptical contour
    3336} pmSourceRadialProfile;
    3437
    3538typedef struct {
    36   float flux;
    37   float fluxErr;
    38   float radius;
    39   float radiusErr;
     39    float flux;
     40    float fluxErr;
     41    float radius;
     42    float radiusErr;
    4043} pmSourceExtendedFlux;
    4144
    4245typedef struct {
    43   pmSourceRadialProfile   *profile;
    44   pmSourceExtendedFlux    *petrosian_50;
    45   pmSourceExtendedFlux    *petrosian_80;
     46    pmSourceRadialFlux     *radFlux;        // raw radial flux information
     47    pmSourceEllipticalFlux *ellipticalFlux; // flux for elliptically-renormalized radii
     48    pmSourceRadialProfile  *radProfile;     // surface brightness profile in specified fixed bins
     49    pmSourceRadialProfile  *petProfile;     // surface brightness profile in petrosian bins
     50    psEllipseAxes axes;                     // shape of elliptical contour
     51    float petrosianFlux;
     52    float petrosianFluxErr;
     53    float petrosianRadius;
     54    float petrosianRadiusErr;
     55    float petrosianR90;
     56    float petrosianR90Err;
     57    float petrosianR50;
     58    float petrosianR50Err;
    4659} pmSourceExtendedPars;
     60
     61pmSourceRadialFlux *pmSourceRadialFluxAlloc();
     62bool psMemCheckSourceRadialFlux(psPtr ptr);
     63
     64pmSourceEllipticalFlux *pmSourceEllipticalFluxAlloc();
     65bool psMemCheckSourceEllipticalFlux(psPtr ptr);
    4766
    4867// *** pmSourceRadialProfile describes the radial profile of a source in elliptical contours, and
    … …  
    5069pmSourceRadialProfile *pmSourceRadialProfileAlloc();
    5170bool psMemCheckSourceRadialProfile(psPtr ptr);
    52 
    53 // *** pmSourceRadialProfileFreeVectors frees the intermediate data values
    54 bool pmSourceRadialProfileFreeVectors(pmSourceRadialProfile *profile);
    5571
    5672// *** pmSourceExtendedPars describes the possible collection of extended flux measurements for a source
    … …  
    6581bool pmSourceRadialProfileSortPair(psVector *index, psVector *extra);
    6682
     83
     84
    6785/// @}
    6886# endif /* PM_SOURCE_H */

  • branches/pap/psModules/src/objects/pmSourceFitSet.c

    r25754 r28003  
    6262    // before pmSourceFitSetInit is called
    6363    for (int i = 0; i < fitSets->n; i++) {
    64         psAssert (fitSets->data[i] == NULL, "failure to init or clear fitSets?");
     64        psAssert(fitSets->data[i] == NULL, "failure to init or clear fitSets?");
    6565    }
    6666    return true;
    … …  
    6868
    6969void pmSourceFitSetDone (void) {
    70     psFree (fitSets);
     70    psFree(fitSets);
     71    fitSets = NULL;
    7172}
    7273
    … …  
    113114}
    114115
    115 pmSourceFitSetData *pmSourceFitSetDataSet (psArray *modelSet) {
    116 
    117     pmSourceFitSetData *thisSet = NULL;
    118 
    119     psAssert (fitSets, "pmSourceFitSetInit not called");
     116pmSourceFitSetData *pmSourceFitSetDataSet (psArray *modelSet)
     117{
     118    psAssert(fitSets, "pmSourceFitSetInit not called");
    120119
    121120    // find the fitSet used by this thread
    122121    pthread_t id = pthread_self();
    123122
    124     // is our ID is already on the stack, abort.
    125     // we do not need to lock to do this....
     123    // If our ID is already on the stack, abort.
     124    // We do need to lock on this because someone might pull one of the fitSets out from under us
     125    pthread_mutex_lock(&fitSetInitMutex);
    126126    for (int i = 0; i < fitSets->n; i++) {
    127         thisSet = fitSets->data[i];
     127        pmSourceFitSetData *thisSet = fitSets->data[i];
    128128        if (!thisSet) continue;
    129129        if (thisSet->thread == id) {
    130             break;
    131         }
    132         thisSet = NULL;
    133     }
    134     psAssert (thisSet == NULL, "pmSourceFitSetDataSet() called but previous entry not cleared");
    135 
    136     thisSet = pmSourceFitSetDataAlloc(modelSet);
    137 
    138     pthread_mutex_lock (&fitSetInitMutex);
    139 
    140     // find an entry that is NULL and place it there
     130            psAbort("pmSourceFitSetDataSet() called but previous entry not cleared");
     131        }
     132    }
     133
     134    // Find an open slot
    141135    for (int i = 0; i < fitSets->n; i++) {
    142136        if (fitSets->data[i]) continue;
    143         fitSets->data[i] = thisSet;
    144         pthread_mutex_unlock (&fitSetInitMutex);
     137        pmSourceFitSetData *thisSet = fitSets->data[i] = pmSourceFitSetDataAlloc(modelSet);
     138        pthread_mutex_unlock(&fitSetInitMutex);
    145139        return thisSet;
    146140    }
    147     pthread_mutex_unlock (&fitSetInitMutex);
    148     psAbort ("no empty slot for new pmSourceFitSetData");
     141    pthread_mutex_unlock(&fitSetInitMutex);
     142    psAbort("no empty slot for new pmSourceFitSetData");
    149143    return NULL;
    150144}
    … …  
    152146pmSourceFitSetData *pmSourceFitSetDataGet (void) {
    153147
    154     psAssert (fitSets, "pmSourceFitSetInit not called");
    155 
    156     // find the fitSet used by this thread
     148    psAssert(fitSets, "pmSourceFitSetInit not called");
     149
     150    // Find the fitSet used by this thread
     151    // We do need to lock on this because someone might pull one of the fitSets out from under us
    157152    pthread_t id = pthread_self();
    158 
    159     // can we find our fitSet?  we do not need to lock to do this....
    160     pmSourceFitSetData *thisSet = NULL;
     153    pthread_mutex_lock(&fitSetInitMutex);
    161154    for (int i = 0; i < fitSets->n; i++) {
    162         thisSet = fitSets->data[i];
     155        pmSourceFitSetData *thisSet = fitSets->data[i];
    163156        if (!thisSet) continue;
    164157        if (thisSet->thread == id) {
    165             break;
    166         }
    167         thisSet = NULL;
    168     }
    169     psAssert (thisSet != NULL, "pmSourceFitSetDataGet() called, but no entry found");
    170 
    171     return thisSet;
    172 }
    173 
    174 void pmSourceFitSetDataClear (void) {
    175 
    176     int i;
    177 
     158            pthread_mutex_unlock(&fitSetInitMutex);
     159            return thisSet;
     160        }
     161    }
     162    psAbort("pmSourceFitSetDataGet() called, but no entry found");
     163}
     164
     165void pmSourceFitSetDataClear (void)
     166{
    178167    psAssert (fitSets, "pmSourceFitSetInit not called");
    179168
    180     // find the fitSet used by this thread
     169    // Find the fitSet used by this thread
     170    // We do need to lock on this because someone might pull one of the fitSets out from under us
    181171    pthread_t id = pthread_self();
    182 
    183     // can we find our fitSet?  we do not need to lock to do this....
    184     pmSourceFitSetData *thisSet = NULL;
    185     for (i = 0; i < fitSets->n; i++) {
    186         thisSet = fitSets->data[i];
     172    pthread_mutex_lock(&fitSetInitMutex);
     173    for (int i = 0; i < fitSets->n; i++) {
     174        pmSourceFitSetData *thisSet = fitSets->data[i];
    187175        if (!thisSet) continue;
    188176        if (thisSet->thread == id) {
    189             break;
    190         }
    191         thisSet = NULL;
    192     }
    193     psAssert (thisSet != NULL, "pmSourceFitSetDataClear() called, but no entry found");
    194 
    195     psFree (thisSet);
    196     fitSets->data[i] = NULL;
     177            psFree(thisSet);
     178            fitSets->data[i] = NULL;
     179            pthread_mutex_unlock(&fitSetInitMutex);
     180            return;
     181        }
     182    }
     183    psAbort("pmSourceFitSetDataClear() called, but no entry found");
     184
    197185    return;
    198186}
    … …  
    481469
    482470            // Convert i/j to image space:
    483             // 0.5 PIX: the coordinate values must be in pixel coords, not index           
     471            // 0.5 PIX: the coordinate values must be in pixel coords, not index
    484472            coord->data.F32[0] = (psF32) (j + 0.5 + source->pixels->col0);
    485473            coord->data.F32[1] = (psF32) (i + 0.5 + source->pixels->row0);

  • branches/pap/psModules/src/objects/pmSourceIO.c

    r27531 r28003  
    551551                }
    552552                if (!strcmp (exttype, "PS1_V1")) {
    553                     status &= pmSourcesWrite_CMF_PS1_V1_XSRC (file->fits, sources, xsrcname, recipe);
     553                    status &= pmSourcesWrite_CMF_PS1_V1_XSRC (file->fits, readout, sources, file->header, xsrcname, recipe);
    554554                }
    555555                if (!strcmp (exttype, "PS1_V2")) {
    556                     status &= pmSourcesWrite_CMF_PS1_V2_XSRC (file->fits, sources, xsrcname, recipe);
     556                    status &= pmSourcesWrite_CMF_PS1_V2_XSRC (file->fits, readout, sources, file->header, xsrcname, recipe);
    557557                }
    558558                if (!strcmp (exttype, "PS1_DV1")) {
    559                     status &= pmSourcesWrite_CMF_PS1_DV1_XSRC (file->fits, sources, xsrcname, recipe);
     559                    status &= pmSourcesWrite_CMF_PS1_DV1_XSRC (file->fits, readout, sources, file->header, xsrcname, recipe);
    560560                }
    561561            }
    … …  
    568568                }
    569569                if (!strcmp (exttype, "PS1_V1")) {
    570                     status &= pmSourcesWrite_CMF_PS1_V1_XFIT (file->fits, sources, xfitname);
     570                    status &= pmSourcesWrite_CMF_PS1_V1_XFIT (file->fits, readout, sources, xfitname);
    571571                }
    572572                if (!strcmp (exttype, "PS1_V2")) {
    573                     status &= pmSourcesWrite_CMF_PS1_V2_XFIT (file->fits, sources, xfitname);
     573                    status &= pmSourcesWrite_CMF_PS1_V2_XFIT (file->fits, readout, sources, xfitname);
    574574                }
    575575                if (!strcmp (exttype, "PS1_DV1")) {
    576                     status &= pmSourcesWrite_CMF_PS1_DV1_XFIT (file->fits, sources, xfitname);
     576                    status &= pmSourcesWrite_CMF_PS1_DV1_XFIT (file->fits, readout, sources, xfitname);
    577577                }
    578578            }

  • branches/pap/psModules/src/objects/pmSourceIO.h

    r27531 r28003  
    3636
    3737bool pmSourcesWrite_CMF_PS1_V1 (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, psMetadata *tableHeader, char *extname);
    38 bool pmSourcesWrite_CMF_PS1_V1_XSRC (psFits *fits, psArray *sources, char *extname, psMetadata *recipe);
    39 bool pmSourcesWrite_CMF_PS1_V1_XFIT (psFits *fits, psArray *sources, char *extname);
     38bool pmSourcesWrite_CMF_PS1_V1_XSRC (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname, psMetadata *recipe);
     39bool pmSourcesWrite_CMF_PS1_V1_XFIT (psFits *fits, pmReadout *readout, psArray *sources, char *extname);
    4040
    4141bool pmSourcesWrite_CMF_PS1_V2 (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, psMetadata *tableHeader, char *extname);
    42 bool pmSourcesWrite_CMF_PS1_V2_XSRC (psFits *fits, psArray *sources, char *extname, psMetadata *recipe);
    43 bool pmSourcesWrite_CMF_PS1_V2_XFIT (psFits *fits, psArray *sources, char *extname);
     42bool pmSourcesWrite_CMF_PS1_V2_XSRC (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname, psMetadata *recipe);
     43bool pmSourcesWrite_CMF_PS1_V2_XFIT (psFits *fits, pmReadout *readout, psArray *sources, char *extname);
    4444
    4545bool pmSourcesWrite_CMF_PS1_DV1 (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, psMetadata *tableHeader, char *extname);
    46 bool pmSourcesWrite_CMF_PS1_DV1_XSRC (psFits *fits, psArray *sources, char *extname, psMetadata *recipe);
    47 bool pmSourcesWrite_CMF_PS1_DV1_XFIT (psFits *fits, psArray *sources, char *extname);
     46bool pmSourcesWrite_CMF_PS1_DV1_XSRC (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname, psMetadata *recipe);
     47bool pmSourcesWrite_CMF_PS1_DV1_XFIT (psFits *fits, pmReadout *readout, psArray *sources, char *extname);
    4848
    4949bool pmSource_CMF_WritePHU (const pmFPAview *view, pmFPAfile *file, pmConfig *config);

  • branches/pap/psModules/src/objects/pmSourceIO_CMF_PS1_DV1.c

    r27531 r28003  
    364364
    365365// XXX this layout is still the same as PS1_DEV_1
    366 bool pmSourcesWrite_CMF_PS1_DV1_XSRC (psFits *fits, psArray *sources, char *extname, psMetadata *recipe)
     366bool pmSourcesWrite_CMF_PS1_DV1_XSRC (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname, psMetadata *recipe)
    367367{
    368368
    … …  
    536536
    537537// XXX this layout is still the same as PS1_DEV_1
    538 bool pmSourcesWrite_CMF_PS1_DV1_XFIT (psFits *fits, psArray *sources, char *extname)
     538bool pmSourcesWrite_CMF_PS1_DV1_XFIT (psFits *fits, pmReadout *readout, psArray *sources, char *extname)
    539539{
    540540

  • branches/pap/psModules/src/objects/pmSourceIO_CMF_PS1_V1.c

    r27177 r28003  
    339339
    340340// XXX this layout is still the same as PS1_DEV_1
    341 bool pmSourcesWrite_CMF_PS1_V1_XSRC (psFits *fits, psArray *sources, char *extname, psMetadata *recipe)
     341bool pmSourcesWrite_CMF_PS1_V1_XSRC (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname, psMetadata *recipe)
    342342{
    343343
    … …  
    516516
    517517// XXX this layout is still the same as PS1_DEV_1
    518 bool pmSourcesWrite_CMF_PS1_V1_XFIT (psFits *fits, psArray *sources, char *extname)
     518bool pmSourcesWrite_CMF_PS1_V1_XFIT (psFits *fits, pmReadout *readout, psArray *sources, char *extname)
    519519{
    520520

  • branches/pap/psModules/src/objects/pmSourceIO_CMF_PS1_V2.c

    r27177 r28003  
    4545// followed by a zero-size matrix, followed by the table data
    4646
    47 bool pmSourcesWrite_CMF_PS1_V2 (psFits *fits, pmReadout *readout, psArray *sources,
    48                                 psMetadata *imageHeader, psMetadata *tableHeader, char *extname)
     47bool pmSourcesWrite_CMF_PS1_V2 (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, psMetadata *tableHeader, char *extname)
    4948{
    5049    PS_ASSERT_PTR_NON_NULL(fits, false);
    … …  
    5453    psArray *table;
    5554    psMetadata *row;
    56     int i;
    5755    psF32 *PAR, *dPAR;
    5856    psEllipseAxes axes;
    … …  
    8381        pmSource *source = (pmSource *) sources->data[0];
    8482        if (source->seq == -1) {
    85           // let's write these out in S/N order
    86           sources = psArraySort (sources, pmSourceSortBySN);
     83            // let's write these out in S/N order
     84            sources = psArraySort (sources, pmSourceSortBySN);
    8785        } else {
    88           sources = psArraySort (sources, pmSourceSortBySeq);
     86            sources = psArraySort (sources, pmSourceSortBySeq);
    8987        }
    9088    }
    … …  
    9492    // we write out PSF-fits for all sources, regardless of quality.  the source flags tell us the state
    9593    // by the time we call this function, all values should be assigned.  let's use asserts to be sure in some cases.
    96     for (i = 0; i < sources->n; i++) {
     94    for (int i = 0; i < sources->n; i++) {
    9795        pmSource *source = (pmSource *) sources->data[i];
    9896
    … …  
    102100        // generated on Alloc, and would thus be wrong for read in sources.
    103101        if (source->seq == -1) {
    104           source->seq = i;
     102            source->seq = i;
    105103        }
    106104
    … …  
    114112            yPos = PAR[PM_PAR_YPOS];
    115113            if (source->mode & PM_SOURCE_MODE_NONLINEAR_FIT) {
    116               xErr = dPAR[PM_PAR_XPOS];
    117               yErr = dPAR[PM_PAR_YPOS];
     114                xErr = dPAR[PM_PAR_XPOS];
     115                yErr = dPAR[PM_PAR_YPOS];
    118116            } else {
    119               // in linear-fit mode, there is no error on the centroid
    120               xErr = source->peak->dx;
    121               yErr = source->peak->dy;
     117                // in linear-fit mode, there is no error on the centroid
     118                xErr = source->peak->dx;
     119                yErr = source->peak->dy;
    122120            }
    123121            if (isfinite(PAR[PM_PAR_SXX]) && isfinite(PAR[PM_PAR_SXX]) && isfinite(PAR[PM_PAR_SXX])) {
    … …  
    212210    }
    213211
     212    // XXX why do we make a copy here to be supplemented with the masks?  why not do this in the calling function?
    214213    psMetadata *header = psMetadataCopy(NULL, tableHeader);
    215214    pmSourceMasksHeader(header);
    … …  
    344343}
    345344
    346 // XXX this layout is still the same as PS1_DEV_1
    347 bool pmSourcesWrite_CMF_PS1_V2_XSRC (psFits *fits, psArray *sources, char *extname, psMetadata *recipe)
     345bool pmSourcesWrite_CMF_PS1_V2_XSRC (psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname, psMetadata *recipe)
    348346{
    349347
    … …  
    354352    psF32 xPos, yPos;
    355353    psF32 xErr, yErr;
     354    int nRow = -1;
    356355
    357356    // create a header to hold the output data
    … …  
    361360    psMetadataAddStr (outhead, PS_LIST_TAIL, "EXTNAME", PS_META_REPLACE, "xsrc table extension", extname);
    362361
     362    pmChip *chip = readout->parent->parent;
     363    pmFPA  *fpa  = chip->parent;
     364
     365    // zero point corrections
     366    bool status1 = false;
     367    bool status2 = false;
     368    float magOffset = 0.0;
     369    float exptime   = psMetadataLookupF32(&status1, fpa->concepts, "FPA.EXPOSURE");
     370    float zeropt    = psMetadataLookupF32(&status2, fpa->concepts, "FPA.ZP");
     371    if (!isfinite(zeropt)) {
     372        zeropt    = psMetadataLookupF32 (&status2, imageHeader, "ZPT_OBS");
     373    }
     374    if (status1 && status2 && (exptime > 0.0)) {
     375        magOffset = zeropt + 2.5*log10(exptime);
     376    }
     377
    363378    // let's write these out in S/N order
    364379    sources = psArraySort (sources, pmSourceSortBySN);
    … …  
    367382
    368383    // which extended source analyses should we perform?
    369     // bool doPetrosian    = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_PETROSIAN");
     384    bool doAnnuli       = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ANNULI");
     385    bool doPetrosian    = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_PETROSIAN");
    370386    // bool doIsophotal    = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ISOPHOTAL");
    371     // bool doAnnuli       = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ANNULI");
    372387    // bool doKron         = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_KRON");
    373388
    374     psVector *radialBinsLower = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.LOWER");
    375     psVector *radialBinsUpper = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.UPPER");
    376     assert (radialBinsLower->n == radialBinsUpper->n);
     389    psVector *radMin = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.LOWER");
     390    psVector *radMax = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.UPPER");
     391    psAssert (radMin->n == radMax->n, "inconsistent annular bins");
     392
     393    // int nRadialBins = 0;
     394    // if (doAnnuli) {
     395    //  // get the max count of radial bins
     396    //  for (int i = 0; i < sources->n; i++) {
     397    //      pmSource *source = sources->data[i];
     398    //      if (!source->extpars) continue;
     399    //      if (!source->extpars->radProfile ) continue;
     400    //         if (!source->extpars->radProfile->binSB) continue;
     401    //      nRadialBins = PS_MAX(nRadialBins, source->extpars->radProfile->binSB->n);
     402    //  }
     403    // }
    377404
    378405    // we write out all sources, regardless of quality.  the source flags tell us the state
    … …  
    410437        psMetadataAdd (row, PS_LIST_TAIL, "Y_EXT_SIG",        PS_DATA_F32, "Sigma in EXT y coordinate",                  yErr);
    411438
    412 # if (0)
     439        float AxialRatio = NAN;
     440        float AxialTheta = NAN;
     441        pmSourceExtendedPars *extpars = source->extpars;
     442        if (extpars) {
     443            AxialRatio = extpars->axes.minor / extpars->axes.major;
     444            AxialTheta = extpars->axes.theta;
     445        }
     446        psMetadataAdd (row, PS_LIST_TAIL, "F25_ARATIO",       PS_DATA_F32, "Axial Ratio of radial profile",              AxialRatio);
     447        psMetadataAdd (row, PS_LIST_TAIL, "F25_THETA",        PS_DATA_F32, "Angle of radial profile ellipse",                  AxialTheta);
     448
    413449        // Petrosian measurements
    414450        // XXX insert header data: petrosian ref radius, flux ratio
     451        // XXX check flags to see if Pet was measured
    415452        if (doPetrosian) {
    416             pmSourcePetrosianValues *petrosian = source->extpars->petrosian;
    417             if (petrosian) {
    418                 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG",        PS_DATA_F32, "Petrosian Magnitude",       petrosian->mag);
    419                 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG_ERR",    PS_DATA_F32, "Petrosian Magnitude Error", petrosian->magErr);
    420                 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS",     PS_DATA_F32, "Petrosian Radius",          petrosian->rad);
    421                 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_ERR", PS_DATA_F32, "Petrosian Radius Error",    petrosian->radErr);
     453            pmSourceExtendedPars *extpars = source->extpars;
     454            if (extpars) {
     455                // XXX note that this mag is either calibrated or instrumental depending on existence of zero point
     456                float mag = (extpars->petrosianFlux > 0.0) ? -2.5*log10(extpars->petrosianFlux) + magOffset : NAN; // XXX zero point
     457                float magErr = (extpars->petrosianFlux > 0.0) ? extpars->petrosianFlux / extpars->petrosianFluxErr : NAN; // XXX zero point
     458                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG",        PS_DATA_F32, "Petrosian Magnitude", mag);
     459                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG_ERR",    PS_DATA_F32, "Petrosian Magnitude Error", magErr);
     460                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS",     PS_DATA_F32, "Petrosian Radius (pix)", extpars->petrosianRadius);
     461                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_ERR", PS_DATA_F32, "Petrosian Radius Error (pix)", extpars->petrosianRadiusErr);
     462                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_50",     PS_DATA_F32, "Petrosian R50 (pix)", extpars->petrosianR50);
     463                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_50_ERR", PS_DATA_F32, "Petrosian R50 Error (pix)", extpars->petrosianR50Err);
     464                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90",     PS_DATA_F32, "Petrosian R90 (pix)", extpars->petrosianR90);
     465                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90_ERR", PS_DATA_F32, "Petrosian R90 Error (pix)", extpars->petrosianR90Err);
    422466            } else {
    423467                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG",        PS_DATA_F32, "Petrosian Magnitude",       NAN);
    … …  
    425469                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS",     PS_DATA_F32, "Petrosian Radius",          NAN);
    426470                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_ERR", PS_DATA_F32, "Petrosian Radius Error",    NAN);
     471                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_50",     PS_DATA_F32, "Petrosian R50 (pix)", NAN);
     472                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_50_ERR", PS_DATA_F32, "Petrosian R50 Error (pix)",NAN);
     473                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90",     PS_DATA_F32, "Petrosian R90 (pix)", NAN);
     474                psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90_ERR", PS_DATA_F32, "Petrosian R90 Error (pix)",NAN);
    427475            }
    428476        }
    429477
     478# if (0)
    430479        // Kron measurements
    431480        if (doKron) {
    … …  
    460509            }
    461510        }
    462 
    463         // Flux Annuli
     511# endif
     512
     513        // Flux Annuli (if we have extended source measurements, we have these.  only optionally save them)
    464514        if (doAnnuli) {
    465             pmSourceAnnuli *annuli = source->extpars->annuli;
    466             if (annuli) {
    467                 psVector *fluxVal = annuli->flux;
    468                 psVector *fluxErr = annuli->fluxErr;
    469                 psVector *fluxVar = annuli->fluxVar;
    470 
    471                 for (int j = 0; j < fluxVal->n; j++) {
    472                     char name[32];
    473                     sprintf (name, "FLUX_VAL_R_%02d", j);
    474                     psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux value in annulus", fluxVal->data.F32[j]);
    475                     sprintf (name, "FLUX_ERR_R_%02d", j);
    476                     psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux error in annulus", fluxErr->data.F32[j]);
    477                     sprintf (name, "FLUX_VAR_R_%02d", j);
    478                     psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux stdev in annulus", fluxVar->data.F32[j]);
    479                 }
    480             } else {
    481                 for (int j = 0; j < radialBinsLower->n; j++) {
    482                     char name[32];
    483                     sprintf (name, "FLUX_VAL_R_%02d", j);
    484                     psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux value in annulus", NAN);
    485                     sprintf (name, "FLUX_ERR_R_%02d", j);
    486                     psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux error in annulus", NAN);
    487                     sprintf (name, "FLUX_VAR_R_%02d", j);
    488                     psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "flux stdev in annulus", NAN);
    489                 }
    490             }
    491         }
    492 
    493 # endif
    494         psArrayAdd (table, 100, row);
    495         psFree (row);
    496     }
    497 
     515            psVector *radSB   = psVectorAlloc(radMin->n, PS_TYPE_F32);
     516            psVector *radFlux = psVectorAlloc(radMin->n, PS_TYPE_F32);
     517            psVector *radFill = psVectorAlloc(radMin->n, PS_TYPE_F32);
     518            psVectorInit (radSB, NAN);
     519            psVectorInit (radFlux, NAN);
     520            psVectorInit (radFill, NAN);
     521            if (!source->extpars) goto empty_annuli;
     522            if (!source->extpars->radProfile) goto empty_annuli;
     523            if (!source->extpars->radProfile->binSB) goto empty_annuli;
     524            psAssert (source->extpars->radProfile->binSum, "programming error");
     525            psAssert (source->extpars->radProfile->binFill, "programming error");
     526            psAssert (source->extpars->radProfile->binSB->n <= radFlux->n, "inconsistent vector lengths");
     527            psAssert (source->extpars->radProfile->binSum->n <= radFlux->n, "inconsistent vector lengths");
     528            psAssert (source->extpars->radProfile->binFill->n <= radFlux->n, "inconsistent vector lengths");
     529
     530            // copy the data from fluxVal (which is not guaranteed to be the full length) to radFlux
     531            for (int j = 0; j < source->extpars->radProfile->binSB->n; j++) {
     532                radSB->data.F32[j]   = source->extpars->radProfile->binSB->data.F32[j];
     533                radFlux->data.F32[j] = source->extpars->radProfile->binSum->data.F32[j];
     534                radFill->data.F32[j] = source->extpars->radProfile->binFill->data.F32[j];
     535            }
     536
     537        empty_annuli:
     538            psMetadataAdd (row, PS_LIST_TAIL, "PROF_SB", PS_DATA_VECTOR, "mean surface brightness annuli", radSB);
     539            psMetadataAdd (row, PS_LIST_TAIL, "PROF_FLUX", PS_DATA_VECTOR, "flux within annuli", radFlux);
     540            psMetadataAdd (row, PS_LIST_TAIL, "PROF_FILL", PS_DATA_VECTOR, "fill factor of annuli", radFill);
     541            psFree (radSB);
     542            psFree (radFlux);
     543            psFree (radFill);
     544        }
     545        if (nRow < 0) {
     546            nRow = row->list->n;
     547        } else {
     548            psAssert (nRow == row->list->n, "inconsistent row lengths");
     549        }
     550        psArrayAdd (table, 100, row);
     551        psFree (row);
     552    }
     553   
    498554    if (table->n == 0) {
    499         if (!psFitsWriteBlank (fits, outhead, extname)) {
    500             psError(psErrorCodeLast(), false, "Unable to write empty sources file.");
    501             psFree(outhead);
    502             psFree(table);
    503             return false;
    504         }
    505         psFree (outhead);
    506         psFree (table);
    507         return true;
    508     }
    509 
     555        if (!psFitsWriteBlank (fits, outhead, extname)) {
     556            psError(psErrorCodeLast(), false, "Unable to write empty sources file.");
     557            psFree(outhead);
     558            psFree(table);
     559            return false;
     560        }
     561        psFree (outhead);
     562        psFree (table);
     563        return true;
     564    }
     565   
    510566    psTrace ("pmFPAfile", 5, "writing ext data %s\n", extname);
    511567    if (!psFitsWriteTable (fits, outhead, table, extname)) {
    512         psError(psErrorCodeLast(), false, "writing ext data %s\n", extname);
    513         psFree (outhead);
    514         psFree(table);
    515         return false;
     568        psError(psErrorCodeLast(), false, "writing ext data %s\n", extname);
     569        psFree (outhead);
     570    psFree(table);
     571    return false;
    516572    }
    517573    psFree (outhead);
    518574    psFree (table);
    519 
     575   
    520576    return true;
    521577}
    522578
    523579// XXX this layout is still the same as PS1_DEV_1
    524 bool pmSourcesWrite_CMF_PS1_V2_XFIT (psFits *fits, psArray *sources, char *extname)
     580bool pmSourcesWrite_CMF_PS1_V2_XFIT (psFits *fits, pmReadout *readout, psArray *sources, char *extname)
    525581{
    526582

  • branches/pap/psModules/src/psmodules.h

    r27672 r28003  
    1010#include <pmKapaPlots.h>
    1111#include <pmVisual.h>
     12#include <ippStages.h>
     13#include <ippDiffMode.h>
    1214
    1315// XXX the following headers define constructs needed by the elements below
Note: See TracChangeset for help on using the changeset viewer.