Index: /trunk/psastro/src/psastro.c =================================================================== --- /trunk/psastro/src/psastro.c (revision 12491) +++ /trunk/psastro/src/psastro.c (revision 12492) @@ -7,4 +7,7 @@ int main (int argc, char **argv) { + + pmConfig *config = NULL; + psArray *refs = NULL; psTimerStart ("complete"); @@ -17,5 +20,5 @@ // load configuration information - pmConfig *config = psastroArguments (argc, argv); + config = psastroArguments (argc, argv); if (!config) usage (); @@ -41,5 +44,5 @@ // load the reference stars overlapping the data stars - psArray *refs = psastroLoadRefstars (config); + refs = psastroLoadRefstars (config); if (!refs) { psErrorStackPrint(stderr, "failed to load reference data\n"); Index: /trunk/psastro/src/psastro.h =================================================================== --- /trunk/psastro/src/psastro.h (revision 12491) +++ /trunk/psastro/src/psastro.h (revision 12492) @@ -78,4 +78,5 @@ bool psastroDumpRefstars (psArray *refstars, char *filename); bool psastroDumpMatches (pmFPA *fpa, char *filename); +bool psastroDumpStars (psArray *stars, char *filename); bool psastroMosaicSetAstrom_tmp (pmFPA *fpa); Index: /trunk/psastro/src/psastroArguments.c =================================================================== --- /trunk/psastro/src/psastroArguments.c (revision 12491) +++ /trunk/psastro/src/psastroArguments.c (revision 12492) @@ -32,5 +32,5 @@ if ((N = psArgumentGet (argc, argv, "-chip"))) { psArgumentRemove (N, &argc, argv); - psMetadataAddStr (config->arguments, PS_LIST_TAIL, "CHIP_SELECTIONS", PS_META_REPLACE, "", psStringCopy(argv[N])); + psMetadataAddStr (config->arguments, PS_LIST_TAIL, "CHIP_SELECTIONS", PS_META_REPLACE, "", argv[N]); psArgumentRemove (N, &argc, argv); } Index: /trunk/psastro/src/psastroAstromGuess.c =================================================================== --- /trunk/psastro/src/psastroAstromGuess.c (revision 12491) +++ /trunk/psastro/src/psastroAstromGuess.c (revision 12492) @@ -40,5 +40,5 @@ double pixelScale = psMetadataLookupF32 (&status, recipe, "PSASTRO.PIXEL.SCALE"); if (!status) { - psError(PS_ERR_IO, false, "Failed to lookup pixel scale"); + psError(PS_ERR_IO, true, "Failed to lookup pixel scale"); return false; } Index: /trunk/psastro/src/psastroChipAstrom.c =================================================================== --- /trunk/psastro/src/psastroChipAstrom.c (revision 12491) +++ /trunk/psastro/src/psastroChipAstrom.c (revision 12492) @@ -52,4 +52,17 @@ } + char *filename = NULL; + char *chipname = psMetadataLookupStr (NULL, chip->concepts, "CHIP.NAME"); + + psStringAppend (&filename, "rawstars.ch.%s.dat", chipname); + psastroDumpStars (rawstars, filename); + psFree (filename); + filename = NULL; + + psStringAppend (&filename, "refstars.ch.%s.dat", chipname); + psastroDumpStars (refstars, filename); + psFree (filename); + filename = NULL; + // save WCS and analysis metadata in update header psMetadata *updates = psMetadataAlloc(); Index: /trunk/psastro/src/psastroCleanup.c =================================================================== --- /trunk/psastro/src/psastroCleanup.c (revision 12491) +++ /trunk/psastro/src/psastroCleanup.c (revision 12492) @@ -12,6 +12,6 @@ pmConceptsDone (); pmConfigDone (); - // fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, stdout, false), "psastro"); - fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, NULL, false), "psastro"); + fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, stdout, false), "psastro"); + // fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, NULL, false), "psastro"); return; Index: /trunk/psastro/src/psastroDemoDump.c =================================================================== --- /trunk/psastro/src/psastroDemoDump.c (revision 12491) +++ /trunk/psastro/src/psastroDemoDump.c (revision 12492) @@ -1,9 +1,38 @@ # include "psastro.h" + +// this function is used for test purposes (-trace psastro.dump.psastroAstromGuess 1) +bool psastroDumpStars (psArray *stars, char *filename) { + + FILE *f = fopen (filename, "w"); + + for (int i = 0; i < stars->n; i++) { + pmAstromObj *obj = stars->data[i]; + + // write out the upward projections + fprintf (f, "%d %f %f %f %f %f %f %f %f %f\n", i, + obj->sky->r, obj->sky->d, obj->Mag, + obj->TP->x, obj->TP->y, + obj->FP->x, obj->FP->y, + obj->chip->x, obj->chip->y); + } + fclose (f); + return true; +} // this function is used for test purposes (-trace psastro.dump.psastroAstromGuess 1) bool psastroDumpRawstars (psArray *rawstars, pmFPA *fpa, pmChip *chip) { - FILE *f1 = fopen ("rawstars.dat", "w"); - FILE *f2 = fopen ("rawstars.down.dat", "w"); + char *filename = NULL; + char *chipname = psMetadataLookupStr (NULL, chip->concepts, "CHIP.NAME"); + + psStringAppend (&filename, "rawstars.up.%s.dat", chipname); + FILE *f1 = fopen (filename, "w"); + psFree (filename); + filename = NULL; + + psStringAppend (&filename, "rawstars.dn.%s.dat", chipname); + FILE *f2 = fopen (filename, "w"); + psFree (filename); + filename = NULL; for (int i = 0; i < rawstars->n; i++) { Index: /trunk/psastro/src/psastroLoadRefstars.c =================================================================== --- /trunk/psastro/src/psastroLoadRefstars.c (revision 12491) +++ /trunk/psastro/src/psastroLoadRefstars.c (revision 12492) @@ -54,5 +54,5 @@ // XXX set getstar in config? // psStringAppend (&getstarLine, "/home/kiawe/eugene/psconfig/dev.lin64/bin/getstar %s -D CATMODE mef -maglim %f -region %f %f %f %f -o %s", catformat, MAGmax, RAmin, DECmin, RAmax, DECmax, tempFile); - psStringAppend (&getstarLine, "getstar %s -D CATMODE mef -maglim %f -region %f %f %f %f -o %s", catformat, MAGmax, RAmin, DECmin, RAmax, DECmax, tempFile); + psStringAppend (&getstarLine, "/home/kiawe/eugene/psconfig/dev.lin64/bin/getstar %s -D CATMODE mef -maglim %f -region %f %f %f %f -o %s", catformat, MAGmax, RAmin, DECmin, RAmax, DECmax, tempFile); psTrace ("psastro", 3, "%s\n", getstarLine); Index: /trunk/psastro/src/psastroMosaicAstrom.c =================================================================== --- /trunk/psastro/src/psastroMosaicAstrom.c (revision 12491) +++ /trunk/psastro/src/psastroMosaicAstrom.c (revision 12492) @@ -5,15 +5,24 @@ bool psastroMosaicAstrom (pmConfig *config, psArray *refs) { + bool status; + // select the current recipe psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, PSASTRO_RECIPE); if (!recipe) { - psError(PSASTRO_ERR_CONFIG, true, "Can't find PSASTRO recipe!\n"); + psError(PSASTRO_ERR_CONFIG, false, "Can't find PSASTRO recipe!\n"); return false; } + + // physical pixel scale in microns per pixel + double pixelScale = psMetadataLookupF32 (&status, recipe, "PSASTRO.PIXEL.SCALE"); + if (!status) { + psError(PS_ERR_IO, false, "Failed to lookup pixel scale"); + return false; + } // select the input data sources pmFPAfile *input = psMetadataLookupPtr (NULL, config->files, "PSASTRO.INPUT"); if (!input) { - psError(PSASTRO_ERR_CONFIG, true, "Can't find input data!\n"); + psError(PSASTRO_ERR_CONFIG, false, "Can't find input data!\n"); return false; } @@ -27,5 +36,8 @@ // given the existing per-chip astrometry, determine matches between raw and ref stars // is this needed? yes, if we didn't do SingleChip astrometry first - psastroMosaicSetMatch (fpa, recipe, 0); + if (!psastroMosaicSetMatch (fpa, recipe, 0)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to match raw and ref stars for mosaic"); + return false; + } if (psTraceGetLevel("psastro.dump") > 0) { psastroDumpMatches (fpa, "match.0.dat"); } @@ -33,5 +45,8 @@ // modify the chip->toFPA scaling to match knowledge about pixel scale, // then recalculate raw and ref positions - psastroMosaicCommonScale (fpa, recipe); + if (!psastroMosaicCommonScale (fpa, recipe)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to set a common scale for the chips"); + return false; + } if (psTraceGetLevel("psastro.dump") > 0) { psastroDumpMatches (fpa, "match.1.dat"); } @@ -39,32 +54,73 @@ // apply the new distortion terms up and down // refit the per-chip terms with linear fits only - psastroMosaicGradients (fpa, recipe); + if (!psastroMosaicGradients (fpa, recipe)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to measure mosaic gradients"); + return false; + } if (psTraceGetLevel("psastro.dump") > 0) { psastroDumpMatches (fpa, "match.2.dat"); } - psastroMosaicChipAstrom (fpa, recipe, false); + // measure the astrometry for the chips under the distortion term + if (!psastroMosaicChipAstrom (fpa, recipe, false)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to measure chip astrometry in mosaic mode"); + return false; + } if (psTraceGetLevel("psastro.dump") > 0) { psastroDumpMatches (fpa, "match.3.dat"); } - // re-perform the match with a slightly tighter circle - psastroMosaicSetMatch (fpa, recipe, 1); + // do a second pass on the distortion with improved chip positions and rotations + // first, re-perform the match with a slightly tighter circle + if (!psastroMosaicSetMatch (fpa, recipe, 1)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to match raw and ref stars for mosaic (2nd pass)"); + return false; + } + if (!psastroMosaicCommonScale (fpa, recipe)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to set a common scale for the chips (2nd pass)"); + return false; + } + if (!psastroMosaicGradients (fpa, recipe)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to measure mosaic gradients (2nd pass)"); + return false; + } + if (psTraceGetLevel("psastro.dump") > 0) { psastroDumpMatches (fpa, "match.4.dat"); } - // do a second pass on the distortion with improved chip positions and rotations - psastroMosaicCommonScale (fpa, recipe); - psastroMosaicGradients (fpa, recipe); - psastroMosaicChipAstrom (fpa, recipe, false); - psastroMosaicSetMatch (fpa, recipe, 1); - + if (!psastroMosaicChipAstrom (fpa, recipe, false)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to measure chip astrometry in mosaic mode (2nd pass)"); + return false; + } + if (psTraceGetLevel("psastro.dump") > 0) { psastroDumpMatches (fpa, "match.5.dat"); } + // now fit the chips under the common distortion with higher-order terms // first, re-perform the match with a slightly tighter circle - psastroMosaicChipAstrom (fpa, recipe, true); + if (!psastroMosaicSetMatch (fpa, recipe, 1)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to match raw and ref stars for mosaic (3rd pass)"); + return false; + } + if (!psastroMosaicChipAstrom (fpa, recipe, true)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to measure chip astrometry in mosaic mode (3rd pass)"); + return false; + } + if (psTraceGetLevel("psastro.dump") > 0) { psastroDumpMatches (fpa, "match.6.dat"); } - // save WCS and analysis metadata in update header + // save WCS and analysis metadata in update header XXX this is still wrong: the pmFPAExtent + // function returns the FP dimensions in pixel units relative to the 0,0 corner of chip + // 0,0. I need to have a function which loops over all cells, determines the pixel + // dimensions for each cell, converts them to chip pixels, then uses the fpa astrometry + // structures to get the total dimensions of the fpa in fpa units. equivalent to + // pmFPAExtent + psRegion *region = pmFPAExtent (fpa); + region->x0 *= pixelScale; + region->x1 *= pixelScale; + region->y0 *= pixelScale; + region->y1 *= pixelScale; + // XXX also need to add DATE/TIME info and NAXIS1, NAXIS2 psMetadata *updates = psMetadataAlloc(); - psMetadataAddS32 (updates, PS_LIST_TAIL, "NAXIS1", PS_META_REPLACE, "fpa dimensions (um)", 30000); - psMetadataAddS32 (updates, PS_LIST_TAIL, "NAXIS2", PS_META_REPLACE, "fpa dimensions (um)", 30000); - psMetadataAddStr (updates, PS_LIST_TAIL, "DATE-OBS", PS_META_REPLACE, "date", "2003-08-28"); - psMetadataAddStr (updates, PS_LIST_TAIL, "UTC-OBS", PS_META_REPLACE, "date", "7:15:46.47"); + psMetadataAddS32 (updates, PS_LIST_TAIL, "NAXIS1", PS_META_REPLACE, "fpa dimensions (mm)", region->x1 - region->x0); + psMetadataAddS32 (updates, PS_LIST_TAIL, "NAXIS2", PS_META_REPLACE, "fpa dimensions (mm)", region->y1 - region->y0); + psFree (region); - pmAstromWriteBilevelMosaic (updates, fpa, NONLIN_TOL); + if (!pmAstromWriteBilevelMosaic (updates, fpa, NONLIN_TOL)) { + psAbort ("failed to save header terms"); + } + psMetadataAddMetadata (fpa->analysis, PS_LIST_TAIL, "PSASTRO.HEADER", PS_META_REPLACE, "psastro header stats", updates); psFree (updates); Index: /trunk/psastro/src/psastroMosaicGetGrads.c =================================================================== --- /trunk/psastro/src/psastroMosaicGetGrads.c (revision 12491) +++ /trunk/psastro/src/psastroMosaicGetGrads.c (revision 12492) @@ -15,4 +15,6 @@ if (!chip->process || !chip->file_exists) { continue; } + psRegion *region = pmChipExtent (chip); + while ((cell = pmFPAviewNextCell (view, fpa, 1)) != NULL) { psTrace ("psastro", 4, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process); @@ -36,8 +38,9 @@ // measure the local gradients for this set of stars - // XXX update the function prototype to accept an incoming gradient structure to which the new elements are added - grads = pmAstromMeasureGradients (grads, rawstars, refstars, match, recipe); + // the new elements are added to the incoming gradient structure + grads = pmAstromMeasureGradients (grads, rawstars, refstars, match, region, 2, 2); } } + psFree (region); } return (grads); Index: /trunk/psastro/src/psastroMosaicGradients.c =================================================================== --- /trunk/psastro/src/psastroMosaicGradients.c (revision 12491) +++ /trunk/psastro/src/psastroMosaicGradients.c (revision 12492) @@ -16,4 +16,6 @@ if (!chip->process || !chip->file_exists) { continue; } + psRegion *region = pmChipExtent (chip); + while ((cell = pmFPAviewNextCell (view, fpa, 1)) != NULL) { psTrace ("psastro", 4, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process); @@ -37,7 +39,9 @@ // measure the local gradients for this set of stars - gradients = pmAstromMeasureGradients (gradients, rawstars, refstars, match, recipe); + // XXX 2,2 are the number of test boxes on the chip. put this in the recipe + gradients = pmAstromMeasureGradients (gradients, rawstars, refstars, match, region, 2, 2); } } + psFree (region); } @@ -59,7 +63,25 @@ } + // physical pixel scale in microns per pixel (FP is in physical units, chip is in pixels) + double pixelScale = psMetadataLookupF32 (&status, recipe, "PSASTRO.PIXEL.SCALE"); + if (!status) { + psError(PS_ERR_IO, false, "Failed to lookup pixel scale"); + return false; + } + // fit the measured gradients with the telescope distortion model (polynomial order based on toTPA) - pmAstromFitDistortion (fpa, gradients, recipe); - psastroMosaicSetAstrom (fpa); + if (!pmAstromFitDistortion (fpa, gradients, pixelScale)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to fit the distortion terms\n"); + psFree (gradients); + psFree (view); + return false; + } + + if (!psastroMosaicSetAstrom (fpa)) { + psError(PSASTRO_ERR_UNKNOWN, false, "failed to apply mosaic distortion terms\n"); + psFree (gradients); + psFree (view); + return false; + } psFree (gradients); Index: /trunk/psastro/src/psastroMosaicRescaleChips.c =================================================================== --- /trunk/psastro/src/psastroMosaicRescaleChips.c (revision 12491) +++ /trunk/psastro/src/psastroMosaicRescaleChips.c (revision 12492) @@ -1,3 +1,5 @@ # include "psastro.h" + +// XXX this file is no longer used // XXX what is this doing? Index: /trunk/psastro/src/psastroOneChip.c =================================================================== --- /trunk/psastro/src/psastroOneChip.c (revision 12491) +++ /trunk/psastro/src/psastroOneChip.c (revision 12492) @@ -10,4 +10,6 @@ bool status; + pmAstromStats *gridStats = NULL; + pmAstromStats *stats = NULL; // supplied radius is in pixels @@ -25,23 +27,29 @@ REQUIRED_RECIPE_VALUE (int minNstar, "PSASTRO.MIN.NSTAR", S32, "failed to find single-chip min allowed stars\n"); - // find initial offset / rotation - pmAstromStats *gridStats = pmAstromGridMatch (rawstars, refstars, recipe); - if (gridStats == NULL) { - psLogMsg ("psastro", 3, "failed to find a grid match solution\n"); - return false; + // do we need to get a rough initial match? + REQUIRED_RECIPE_VALUE (bool gridSearch, "PSASTRO.GRID.SEARCH", Bool, "failed to find chip grid-search option\n"); + + if (gridSearch) { + // find initial offset / rotation + gridStats = pmAstromGridMatch (rawstars, refstars, recipe); + if (gridStats == NULL) { + psLogMsg ("psastro", 3, "failed to find a grid match solution\n"); + return false; + } + psLogMsg ("psastro", 3, "basic grid search result - offset: %f,%f pixels, rotation: %f deg\n", gridStats->offset.x, gridStats->offset.y, DEG_RAD*gridStats->angle); + + // tweak the position by finding peak of matches stars + stats = pmAstromGridTweak (rawstars, refstars, recipe, gridStats); + if (stats == NULL) { + psLogMsg ("psastro", 3, "failed to measure tweaked grid solution\n"); + psFree (gridStats); + return false; + } + psLogMsg ("psastro", 3, "tweak grid search result - offset: %f,%f pixels, rotation: %f deg\n", stats->offset.x, stats->offset.y, DEG_RAD*stats->angle); + + // adjust the chip.toFPA terms only + pmAstromGridApply (chip->toFPA, stats); + psastroUpdateChipToFPA (fpa, chip, rawstars, refstars); } - psLogMsg ("psastro", 3, "basic grid search result - offset: %f,%f pixels, rotation: %f deg\n", gridStats->offset.x, gridStats->offset.y, DEG_RAD*gridStats->angle); - - // tweak the position by finding peak of matches stars - pmAstromStats *stats = pmAstromGridTweak (rawstars, refstars, recipe, gridStats); - if (stats == NULL) { - psLogMsg ("psastro", 3, "failed to measure tweaked grid solution\n"); - return false; - } - psLogMsg ("psastro", 3, "tweak grid search result - offset: %f,%f pixels, rotation: %f deg\n", stats->offset.x, stats->offset.y, DEG_RAD*stats->angle); - - // adjust the chip.toFPA terms only - pmAstromGridApply (chip->toFPA, stats); - psastroUpdateChipToFPA (fpa, chip, rawstars, refstars); // use small radius to match stars @@ -49,4 +57,6 @@ if (match == NULL) { psLogMsg ("psastro", 3, "failed to find radius-matched sources\n"); + psFree (stats); + psFree (gridStats); return false; } @@ -85,4 +95,8 @@ if (!results) { psLogMsg ("psastro", 3, "failed to perform the matched fit\n"); + psFree (match); + psFree (stats); + psFree (fitStats); + psFree (gridStats); return false; } Index: /trunk/psastro/src/psastroRefstarSubset.c =================================================================== --- /trunk/psastro/src/psastroRefstarSubset.c (revision 12491) +++ /trunk/psastro/src/psastroRefstarSubset.c (revision 12492) @@ -56,5 +56,11 @@ if (psTraceGetLevel("psastro.dump") > 0) { - psastroDumpRefstars (subset, "refstars.subset.dat"); + pmChip *chip = readout->parent->parent; + + char *filename = NULL; + char *chipname = psMetadataLookupStr (NULL, chip->concepts, "CHIP.NAME"); + psStringAppend (&filename, "refstars.%s.dat", chipname); + psastroDumpRefstars (subset, filename); + psFree (filename); } Index: /trunk/psastro/src/psastroUtils.c =================================================================== --- /trunk/psastro/src/psastroUtils.c (revision 12491) +++ /trunk/psastro/src/psastroUtils.c (revision 12492) @@ -79,16 +79,12 @@ bool psastroUpdateChipToFPA (pmFPA *fpa, pmChip *chip, psArray *rawstars, psArray *refstars) { - // XXX this region needs to be defined more intelligently - // psRegion *region = pmChipExtent (chip); - psRegion *region = psRegionAlloc (0, 2000, 0, 2000); - + psRegion *region = pmChipExtent (chip); + region->x1 -= region->x0; + region->y1 -= region->y0; + region->x0 = 0; + region->y0 = 0; psFree (chip->fromFPA); chip->fromFPA = psPlaneTransformInvert (NULL, chip->toFPA, *region, 20); - -// PS_POLY_PRINT_2D(chip->toFPA->x); -// PS_POLY_PRINT_2D(chip->toFPA->y); - -// PS_POLY_PRINT_2D(chip->fromFPA->x); -// PS_POLY_PRINT_2D(chip->fromFPA->y); + psFree (region); for (int i = 0; i < rawstars->n; i++) { @@ -105,5 +101,4 @@ } - psFree (region); return true; }