IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Mar 9, 2009, 3:41:26 PM (17 years ago)
Author:
beaumont
Message:

sync with trunk

Location:
branches/cnb_branches/cnb_branch_20090215
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • branches/cnb_branches/cnb_branch_20090215

  • branches/cnb_branches/cnb_branch_20090215/psModules/src/config/Makefile.am

    r23199 r23239  
    11noinst_LTLIBRARIES = libpsmodulesconfig.la
    22
    3 PSMODULES_VERSION=`if [ -e ../../VERSION ]; then cat ../../VERSION; else svnversion; fi`
    4 PSMODULES_BRANCH=`if [ -e ../../BRANCH ]; then cat ../../BRANCH; else svn info | sed -n -e '/URL:/ h' -e '/Repository Root:/ { x; H; x; s|Repository Root: \(.*\)\nURL: \1\(.*\)|\2| ; s|^/|| ; s|/[a-zA-Z]*/src.*|| ; p }'; fi`
    5 PSMODULES_SOURCE=`if [ -e ../../SOURCE ]; then cat ../../SOURCE; else svn info | sed -n -e 's/Repository UUID: // p'; fi`
     3# PSMODULES_VERSION=`if [ -e ../../VERSION ]; then cat ../../VERSION; else svnversion; fi`
     4# PSMODULES_BRANCH=`if [ -e ../../BRANCH ]; then cat ../../BRANCH; else svn info | sed -n -e '/URL:/ h' -e '/Repository Root:/ { x; H; x; s|Repository Root: \(.*\)\nURL: \1\(.*\)|\2| ; s|^/|| ; s|/[a-zA-Z]*/src.*|| ; p }'; fi`
     5# PSMODULES_SOURCE=`if [ -e ../../SOURCE ]; then cat ../../SOURCE; else svn info | sed -n -e 's/Repository UUID: // p'; fi`
    66
    77# Force recompilation of pmVersion.c, since it gets the version information
    8 pmVersion.c: FORCE
    9         touch pmVersion.c
    10 FORCE: ;
     8# pmVersion.c: FORCE
     9#       touch pmVersion.c
     10# FORCE: ;
    1111
    12 libpsmodulesconfig_la_CPPFLAGS = $(SRCINC) $(PSMODULES_CFLAGS) -DPSMODULES_VERSION=\"$(PSMODULES_VERSION)\" -DPSMODULES_BRANCH=\"$(PSMODULES_BRANCH)\" -DPSMODULES_SOURCE=\"$(PSMODULES_SOURCE)\"
     12libpsmodulesconfig_la_CPPFLAGS = $(SRCINC) $(PSMODULES_CFLAGS) -DPSMODULES_VERSION=$(SVN_VERSION) -DPSMODULES_BRANCH=$(SVN_BRANCH) -DPSMODULES_SOURCE=$(SVN_SOURCE)
    1313libpsmodulesconfig_la_LDFLAGS  = -release $(PACKAGE_VERSION)
    1414libpsmodulesconfig_la_SOURCES  = \

  • branches/cnb_branches/cnb_branch_20090215/psModules/src/config/pmConfig.c

    r21314 r23239  
    11101110// camera.  If we are discovering the camera (config->camera == NULL), then we also load the
    11111111// recipe files for the camera.
    1112 psMetadata *pmConfigCameraFormatFromHeader(psMetadata **camera, psString *cameraName, psString *formatName, 
    1113                                            pmConfig *config, const psMetadata *header, bool readRecipes)
     1112psMetadata *pmConfigCameraFormatFromHeader(psMetadata **camera, psString *cameraName, psString *formatName,
     1113                                           pmConfig *config, const psMetadata *header, bool readRecipes)
    11141114{
    11151115    PS_ASSERT_PTR_NON_NULL(config, NULL);
    … …  
    16721672                // if the object does not exist and create isn't set, then we
    16731673                // should puke
    1674                 psError(PM_ERR_SYS, true, "Unable to access file %s", filename);
     1674                psError(PM_ERR_SYS, true, "Unable to access file %s: %s", filename, nebErr(server));
    16751675                nebServerFree(server);
    16761676                return NULL;
    … …  
    16851685        if (trunc) {
    16861686            if(truncate(path, 0) != 0) {
    1687                 psError(PS_ERR_IO, true, "Failed to truncate Nebulous file %s (real name %s)\n", filename, path);
     1687                psError(PS_ERR_IO, true, "Failed to truncate Nebulous file %s (real name %s)\n",
     1688                        filename, path);
    16881689                return NULL;
    16891690            }

  • branches/cnb_branches/cnb_branch_20090215/psModules/src/config/pmVersion.c

    r23199 r23239  
    88#include "pmVersion.h"
    99
    10 psString psModulesVersion(void)
    11 {
    1210#ifndef PSMODULES_VERSION
    1311#error "PSMODULES_VERSION is not set"
    … …  
    1614#error "PSMODULES_BRANCH is not set"
    1715#endif
    18     return psStringCopy(PSMODULES_BRANCH "@" PSMODULES_VERSION);
     16#ifndef PSMODULES_SOURCE
     17#error "PSMODULES_SOURCE is not set"
     18#endif
     19
     20#define xstr(s) str(s)
     21#define str(s) #s
     22
     23psString psModulesVersion(void)
     24{
     25    char *value = NULL;
     26    psStringAppend(&value, "%s@%s", xstr(PSMODULES_BRANCH), xstr(PSMODULES_VERSION));
     27    return value;
    1928}
    2029
    2130psString psModulesSource(void)
    2231{
    23 #ifndef PSMODULES_SOURCE
    24 #error "PSMODULES_SOURCE is not set"
    25 #endif
    26     return psStringCopy(PSMODULES_SOURCE);
     32    return psStringCopy(xstr(PSMODULES_SOURCE));
    2733}
    2834

  • branches/cnb_branches/cnb_branch_20090215/psModules/src/objects/pmSourceMatch.c

    r23199 r23239  
    152152
    153153
    154 psArray *pmSourceMatchSources(const psArray *sourceArrays, float radius)
     154psArray *pmSourceMatchSources(const psArray *sourceArrays, float radius, bool cullSingles)
    155155{
    156156    PS_ASSERT_ARRAY_NON_NULL(sourceArrays, NULL);
    … …  
    220220
    221221        // Match with the master list
    222 
    223         psTree *tree = psTreePlant(2, SOURCES_MAX_LEAF, xMaster, yMaster); // kd Tree with sources
    224         long numMatch = 0;              // Number of matches
    225 
    226         long size = numMaster + numSources; // New size
    227         xMaster = psVectorRealloc(xMaster, size);
    228         yMaster = psVectorRealloc(yMaster, size);
    229         matches = psArrayRealloc(matches, size);
    230 
    231         psVector *coords = psVectorAlloc(2, PS_TYPE_F32); // Coordinates for tree lookup
    232         for (int j = 0; j < numSources; j++) {
    233             coords->data.F32[0] = xImage->data.F32[j];
    234             coords->data.F32[1] = yImage->data.F32[j];
    235             long index = psTreeNearestWithin(tree, coords, radius); // Match index
    236             if (index >= 0) {
    237                 // Record the match
    238                 pmSourceMatch *match = matches->data[index]; // Match data
    239                 pmSourceMatchAdd(match, magImage->data.F32[j], magErrImage->data.F32[j], i, j);
    240                 numMatch++;
    241             } else {
    242                 // Add to the master list
    243                 pmSourceMatch *match = pmSourceMatchAlloc(numImages); // Match data
    244                 pmSourceMatchAdd(match, magImage->data.F32[j], magErrImage->data.F32[j], i, j);
    245                 xMaster->data.F32[numMaster] = xImage->data.F32[j];
    246                 yMaster->data.F32[numMaster] = yImage->data.F32[j];
    247                 matches->data[numMaster] = match;
    248                 numMaster++;
    249                 xMaster->n = yMaster->n = matches->n = numMaster;
    250             }
    251 
    252         }
    253         psFree(coords);
    254         psFree(tree);
    255 
     222        {
     223            psTree *tree = psTreePlant(2, SOURCES_MAX_LEAF, xMaster, yMaster); // kd Tree with sources
     224            long numMatch = 0;              // Number of matches
     225
     226            long size = numMaster + numSources; // New size
     227            xMaster = psVectorRealloc(xMaster, size);
     228            yMaster = psVectorRealloc(yMaster, size);
     229            matches = psArrayRealloc(matches, size);
     230
     231            psVector *coords = psVectorAlloc(2, PS_TYPE_F32); // Coordinates for tree lookup
     232            for (int j = 0; j < numSources; j++) {
     233                coords->data.F32[0] = xImage->data.F32[j];
     234                coords->data.F32[1] = yImage->data.F32[j];
     235                long index = psTreeNearestWithin(tree, coords, radius); // Match index
     236                if (index >= 0) {
     237                    // Record the match
     238                    pmSourceMatch *match = matches->data[index]; // Match data
     239                    pmSourceMatchAdd(match, magImage->data.F32[j], magErrImage->data.F32[j], i, j);
     240                    numMatch++;
     241                } else {
     242                    // Add to the master list
     243                    pmSourceMatch *match = pmSourceMatchAlloc(numImages); // Match data
     244                    pmSourceMatchAdd(match, magImage->data.F32[j], magErrImage->data.F32[j], i, j);
     245                    xMaster->data.F32[numMaster] = xImage->data.F32[j];
     246                    yMaster->data.F32[numMaster] = yImage->data.F32[j];
     247                    matches->data[numMaster] = match;
     248                    numMaster++;
     249                    xMaster->n = yMaster->n = matches->n = numMaster;
     250                }
     251
     252            }
     253            psFree(coords);
     254            psFree(tree);
     255        }
    256256
    257257    match_image_done:
    … …  
    263263    }
    264264
    265     // Now cull the matches that contain only a single star
    266     int numGood = 0;                    // Number of good matches
     265    if (cullSingles) {
     266        // Now cull the matches that contain only a single star
     267        int numGood = 0;                    // Number of good matches
     268        for (int i = 0; i < matches->n; i++) {
     269            pmSourceMatch *match = matches->data[i]; // Match of interest
     270            if (match->num > 1) {
     271                if (i != numGood) {
     272                    psFree(matches->data[numGood]);
     273                    matches->data[numGood] = psMemIncrRefCounter(match);
     274                }
     275                numGood++;
     276            }
     277        }
     278        matches->n = numGood;
     279        for (int i = numGood; i < numMaster; i++) {
     280            psFree(matches->data[i]);
     281            matches->data[i] = NULL;
     282        }
     283    }
     284
     285    return matches;
     286}
     287
     288
     289psArray *pmSourceMatchMerge(psArray *sourceArrays, float radius)
     290{
     291    PS_ASSERT_ARRAY_NON_NULL(sourceArrays, NULL);
     292    PS_ASSERT_FLOAT_LARGER_THAN(radius, 0.0, NULL);
     293
     294    psArray *matches = pmSourceMatchSources(sourceArrays, radius, false); // Source matches
     295    if (!matches) {
     296        psError(PS_ERR_UNKNOWN, false, "Unable to match source lists.");
     297        return NULL;
     298    }
     299
     300    int index = 0;                      // Index to current position on list
    267301    for (int i = 0; i < matches->n; i++) {
     302        pmSource *source = NULL;        // Source to put in merged list
    268303        pmSourceMatch *match = matches->data[i]; // Match of interest
    269         if (match->num > 1) {
    270             if (i != numGood) {
    271                 psFree(matches->data[numGood]);
    272                 matches->data[numGood] = psMemIncrRefCounter(match);
    273                 //                psFree(match);
    274             }
    275             numGood++;
    276         }
    277     }
    278     matches->n = numGood;
    279     for (int i = numGood; i < numMaster; i++) {
    280         psFree(matches->data[i]);
    281         matches->data[i] = NULL;
    282     }
     304        for (int j = 0; j < match->num && !source; j++) {
     305            if (!isfinite(match->mag->data.F32[j]) || !isfinite(match->magErr->data.F32[j])) {
     306                continue;
     307            }
     308            int image = match->image->data.S32[j]; // Index of image
     309            int index = match->index->data.S32[j]; // Index of source for image
     310            psArray *list = sourceArrays->data[image]; // List of interest
     311            source = list->data[index];
     312            break;
     313        }
     314
     315        if (source) {
     316            psFree(matches->data[index]);
     317            matches->data[index] = psMemIncrRefCounter(source);
     318            index++;
     319        }
     320    }
     321
     322    // Clear out the rest of the list
     323    int num = index;                    // Number of good sources
     324    for (; index < matches->n; index++) {
     325        psFree(matches->data[index]);
     326        matches->data[index] = NULL;
     327    }
     328    matches->n = num;
    283329
    284330    return matches;
    285331}
    286 
    287332
    288333// Iterate on the star magnitudes and image transparencies
    … …  
    524569psVector *pmSourceMatchRelphot(const psArray *matches, // Array of matches
    525570                               const psVector *zp, // Zero points for each image (including airmass term)
    526                                int maxIter, // Maximum number of iterations
    527571                               float tol, // Relative tolerance for convergence
     572                               int iter1, // Number of iterations for pass 1
     573                               float rej1, // Limit on rejection between iterations for pass 1
     574                               float sys1, // Systematic error in measurements for pass 1
     575                               int iter2, // Number of iterations for pass 2
     576                               float rej2, // Limit on rejection between iterations for pass 2
     577                               float sys2, // Systematic error in measurements for pass 2
    528578                               float rejLimit, // Limit on rejection between iterations
    529579                               int transIter, // Clipping iterations for transparency
    530580                               float transClip, // Clipping level for transparency
    531                                float photoLevel, // Level at which we declare image is photometric
    532                                float starClip, // Clipping for stars
    533                                float sysErr // Systematic error in measurements
     581                               float photoLevel // Level at which we declare image is photometric
    534582                               )
    535583{
    … …  
    539587    PS_ASSERT_FLOAT_LARGER_THAN(transClip, 0.0, NULL);
    540588
    541     sysErr *= sysErr;
     589    sys1 *= sys1;
     590    sys2 *= sys2;
    542591
    543592    int numImages = zp->n;              // Number of images
    … …  
    552601
    553602    float chi2 = sourceMatchRelphotIterate(trans, stars, badImage, matches, zp,
    554                                            photo, sysErr); // chi^2 for solution
     603                                           photo, sys1); // chi^2 for solution
    555604    psTrace("psModules.objects", 1, "Initial: chi^2 = %f\n", chi2);
    556605    float lastChi2 = INFINITY;          // chi^2 on last iteration
    557606    float fracRej = INFINITY;        // Fraction of measurements rejected
    558607
    559     // in the first passes, the transparencies are not well deteremined: use high systematic error and the rejection thresholds
    560     for (int i = 0; i < 5; i++) {
     608    // In the first passes, the transparencies are not well deteremined: use high systematic error and
     609    // rejection thresholds
     610    for (int i = 0; i < iter1; i++) {
    561611
    562612        // Identify photometric nights
    … …  
    572622        psTrace("psModules.objects", 3, "Pass 1: Determined %d/%d are photometric", numPhoto, numImages);
    573623
    574         // XXX use 20 sigma rejection and 0.1 mag systematic error (move these to the recipe)
    575         fracRej = sourceMatchRelphotReject(trans, stars, matches, zp, photo, badImage, 20.0, PS_SQR(0.1));
     624        fracRej = sourceMatchRelphotReject(trans, stars, matches, zp, photo, badImage, rej1, sys1);
    576625        psTrace("psModules.objects", 3, "Pass 1: %f%% of measurements rejected", fracRej * 100);
    577626
    578         // XXX use 0.05 mag systematic error (move these to the recipe)
    579         chi2 = sourceMatchRelphotIterate(trans, stars, badImage, matches, zp, photo, PS_SQR(0.1));
     627        chi2 = sourceMatchRelphotIterate(trans, stars, badImage, matches, zp, photo, sys1);
    580628        psTrace("psModules.objects", 1, "Pass 1: iter = %d: chi^2 = %f rejected = %f\n", i, chi2, fracRej);
    581629    }
    582630
    583     for (int i = 0; i < maxIter && (fabsf(lastChi2 - chi2) > tol * chi2 || fracRej > rejLimit); i++) {
     631    for (int i = 0; i < iter2 && (fabsf(lastChi2 - chi2) > tol * chi2 || fracRej > rejLimit); i++) {
    584632        lastChi2 = chi2;
    585633
    … …  
    594642            return NULL;
    595643        }
    596         psTrace("psModules.objects", 3, "Determined %d/%d are photometric", numPhoto, numImages);
    597 
    598         fracRej = sourceMatchRelphotReject(trans, stars, matches, zp, photo, badImage, starClip, sysErr);
    599         psTrace("psModules.objects", 3, "%f%% of measurements rejected", fracRej * 100);
    600 
    601         chi2 = sourceMatchRelphotIterate(trans, stars, badImage, matches, zp, photo, sysErr);
    602         psTrace("psModules.objects", 1, "iter = %d: chi^2 = %f rejected = %f\n", i, chi2, fracRej);
     644        psTrace("psModules.objects", 3, "Pass 2: Determined %d/%d are photometric", numPhoto, numImages);
     645
     646        fracRej = sourceMatchRelphotReject(trans, stars, matches, zp, photo, badImage, rej2, sys2);
     647        psTrace("psModules.objects", 3, "Pass 2: %f%% of measurements rejected", fracRej * 100);
     648
     649        chi2 = sourceMatchRelphotIterate(trans, stars, badImage, matches, zp, photo, sys2);
     650        psTrace("psModules.objects", 1, "Pass 2: iter = %d: chi^2 = %f rejected = %f\n", i, chi2, fracRej);
    603651    }
    604652

  • branches/cnb_branches/cnb_branch_20090215/psModules/src/objects/pmSourceMatch.h

    r20953 r23239  
    4848///
    4949/// Returns an array of psSourceMatch
    50 psArray *pmSourceMatchSources(const psArray *sourceArrays, // Array of arrays of sources on each image
    51                               float radius // Correlation radius
     50psArray *pmSourceMatchSources(const psArray *sourceArrays, ///< Array of arrays of sources on each image
     51                              float radius, ///< Matching radius
     52                              bool cullSingles ///< Cull "matches" with only a single source?
    5253                              );
     54
     55/// Merge two source lists
     56///
     57/// Sources are pulled from the lists into a new list, with no effort made to adjust them.
     58psArray *pmSourceMatchMerge(psArray *sourceArrays, ///< Array of arrays of sources on each image
     59                            float radius ///< Matching radius
     60    );
    5361
    5462/// Perform relative photometry to calibrate images
    5563psVector *pmSourceMatchRelphot(const psArray *matches, // Array of matches
    5664                               const psVector *zp, // Zero points for each image (including airmass term)
    57                                int maxIter, // Maximum number of iterations
    5865                               float tol, // Relative tolerance for convergence
     66                               int iter1, // Number of iterations for pass 1
     67                               float rej1, // Limit on rejection between iterations for pass 1
     68                               float sys1, // Systematic error in measurements for pass 1
     69                               int iter2, // Number of iterations for pass 2
     70                               float rej2, // Limit on rejection between iterations for pass 2
     71                               float sys2, // Systematic error in measurements for pass 2
    5972                               float rejLimit, // Limit on rejection between iterations
    6073                               int transIter, // Clipping iterations for transparency
    6174                               float transClip, // Clipping level for transparency
    62                                float photoLevel, // Level at which we declare image is photometric
    63                                float starClip, // Clipping for stars
    64                                float sysErr // Systematic error in measurements
     75                               float photoLevel // Level at which we declare image is photometric
    6576    );
    6677
Note: See TracChangeset for help on using the changeset viewer.