Index: /trunk/ppBackground/src/ppBackground.c =================================================================== --- /trunk/ppBackground/src/ppBackground.c (revision 36581) +++ /trunk/ppBackground/src/ppBackground.c (revision 36582) @@ -18,17 +18,19 @@ ppBackgroundErrorRegister(); + printf("Initializing data\n"); ppBackgroundData *data = ppBackgroundDataInit(&argc, argv); if (!data) { goto DIE; } - + printf("Reading Arguments\n"); if (!ppBackgroundArguments(data, argc, argv)) { goto DIE; } - + printf("Setting up Camera\n"); if (!ppBackgroundCamera(data)) { goto DIE; } + printf("Looping over data!\n"); if (!ppBackgroundLoop(data)) { goto DIE; Index: /trunk/ppBackground/src/ppBackgroundStack.c =================================================================== --- /trunk/ppBackground/src/ppBackgroundStack.c (revision 36581) +++ /trunk/ppBackground/src/ppBackgroundStack.c (revision 36582) @@ -18,17 +18,20 @@ ppBackgroundErrorRegister(); + printf("Data init\n"); ppBackgroundStackData *data = ppBackgroundStackDataInit(&argc, argv); if (!data) { goto DIE; } - + printf("Read arguments\n"); if (!ppBackgroundStackArguments(data, argc, argv)) { goto DIE; } + printf("Setup camera\n"); if (!ppBackgroundStackCamera(data)) { goto DIE; } + printf("Do Loop\n"); if (!ppBackgroundStackLoop(data)) { goto DIE; Index: /trunk/ppBackground/src/ppBackgroundStack.h =================================================================== --- /trunk/ppBackground/src/ppBackgroundStack.h (revision 36581) +++ /trunk/ppBackground/src/ppBackgroundStack.h (revision 36582) @@ -23,5 +23,6 @@ psF32 dec_min; psF32 dec_max; - + + psArray *stack_data; psArray *stacks; // List of stacks to be corrected. psString outRoot; // Output root name Index: /trunk/ppBackground/src/ppBackgroundStackCamera.c =================================================================== --- /trunk/ppBackground/src/ppBackgroundStackCamera.c (revision 36581) +++ /trunk/ppBackground/src/ppBackgroundStackCamera.c (revision 36582) @@ -36,6 +36,7 @@ psString stackName = data->stacks->data[0]; fileArguments("IMAGE", stackName, "Input image", data->config); - pmFPAfile *stack = pmFPAfileDefineFromArgs(&status, data->config, "PPBACKGROUND.STACK", - "IMAGE"); + psFree(stackName); + pmFPAfile *stack = pmFPAfileDefineFromArgs(&status, data->config, + "PPBACKGROUND.STACK", "IMAGE"); if (!status || !stack) { psError(psErrorCodeLast(), false, "Failed to build file from PPBACKGROUND.STACK"); @@ -43,9 +44,16 @@ } + if (!pmAstromReadBilevelMosaic(stack->fpa,stack->fpa->hdu->header)) { + psError(psErrorCodeLast(), false, "Unable to read WCS astrometry for input stack."); + return false; + } + psArrayAdd(data->stack_data,data->stack_data->n,psMemIncrRefCounter(stack)); + // Read over the input background models. psMetadataIterator *iter = psMetadataIteratorAlloc(data->contents, PS_LIST_HEAD, NULL); // Iterator psMetadataItem *item; // Item from iteration - // int i = 0; + int i = 0; while ((item = psMetadataGetAndIncrement(iter))) { + i++; if (item->type != PS_DATA_METADATA) { psError(PPBACKGROUND_ERR_ARGUMENTS, true, @@ -58,5 +66,5 @@ psString smfFileName = psMetadataLookupStr(NULL, input, "astrom"); - + psMetadata *modelContent = psMetadataLookupMetadata(NULL, input, "models"); // Read the smf file from this item fileArguments("astrom",smfFileName,"",config); @@ -68,5 +76,5 @@ return NULL; } - + printf("CZW: Item %d\n",i); // find the FPA phu bool bilevelAstrometry = false; @@ -91,15 +99,12 @@ } const char *chipName = psMetadataLookupStr(NULL, chip->concepts, "CHIP.NAME"); // Name of chip -/* if (data->chipName && strcmp(chipName, data->chipName) != 0) { */ -/* continue; */ -/* } */ if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { psError(psErrorCodeLast(), false, "Error loading data from files."); return false; } - if (chip->cells->n != 1) { psWarning("More than one cell present for chip %d", view->chip); } + // read WCS data from the corresponding header pmHDU *hdu = pmFPAviewThisHDU (view, smfFile->fpa); @@ -124,13 +129,10 @@ // FIX Load model data for this chip - psString modelFileName = psMetadataLookupStr(NULL, input, chipName); + psString modelFileName = pmConfigConvertFilename(psMetadataLookupStr(NULL, modelContent, chipName), + config, PM_FPA_MODE_READ, false); psFits *modelFits = psFitsOpen(modelFileName,"r"); - int nC, nR; - psElemType *type; - psRegion region; - psFitsImageSize(&nC,&nR,type,modelFits,region); - psImage *image = psFitsReadImage(modelFits,region,0); + psImage *image = psFitsReadImage(modelFits,psRegionSet(0,0,0,0),0); psMetadata *header = psFitsReadHeader(NULL,modelFits); - + // Allocate the data structures for this chip psImage *raim = psImageAlloc(image->numCols,image->numRows,PS_TYPE_F32); @@ -177,18 +179,37 @@ // Attach vectors to teh structure of the chip - psMetadataAddPtr(chip->analysis,PS_LIST_TAIL,"bkg image", PS_DATA_UNKNOWN | PS_META_REPLACE, "ota space X vector", image); - psMetadataAddPtr(chip->analysis,PS_LIST_TAIL,"bkg ra", PS_DATA_UNKNOWN | PS_META_REPLACE, "ota space X vector", raim); - psMetadataAddPtr(chip->analysis,PS_LIST_TAIL,"bkg dec", PS_DATA_UNKNOWN | PS_META_REPLACE, "ota space X vector", decim); - psMetadataAddPtr(chip->analysis,PS_LIST_TAIL,"bkg calibrated data", PS_DATA_UNKNOWN | PS_META_REPLACE, "ota space corrected data", model); + psMetadataAddImage(chip->analysis,PS_LIST_TAIL, + "bkg image", 0, + "ota space X vector", image); + psMetadataAddImage(chip->analysis,PS_LIST_TAIL, + "bkg ra", 0, + "ota space ra vector", raim); + psMetadataAddImage(chip->analysis,PS_LIST_TAIL, + "bkg dec", 0, + "ota space dec vector", decim); + psMetadataAddImage(chip->analysis,PS_LIST_TAIL, + "bkg calibrated data", 0, + "ota space corrected data", model); // Define default background model parameters, using the assumption: // observed = camera + offset + scale * astrophysical - psMetadataAddF32(chip->analysis,PS_LIST_TAIL,"bkg offset", PS_META_REPLACE, "background offset for this exposure/ota pair", 0.0); - psMetadataAddF32(chip->analysis,PS_LIST_TAIL,"bkg scale", PS_META_REPLACE, "background scale parameter for this exposure/ota pair", 1.0); + psMetadataAddF32(chip->analysis,PS_LIST_TAIL, + "bkg offset", PS_META_REPLACE, + "background offset for this exposure/ota pair", 0.0); + psMetadataAddF32(chip->analysis,PS_LIST_TAIL, + "bkg scale", PS_META_REPLACE, + "background scale parameter for this exposure/ota pair", 1.0); // FIX Free model data for this chip + psFree(modelFileName); + psFree(modelFits); psFree(header); - psFree(modelFits); - psFree(modelFileName); - + psFree(pix); + psFree(fp); + psFree(tp); + psFree(sky); + psFree(image); + psFree(raim); + psFree(decim); + psFree(model); // Check to see if we've loaded or allocated an OTA solution container for this chip @@ -197,5 +218,7 @@ if (data->fit_OTAS) { // We are fitting OTAs, so allocate a new image psImage *solution = psImageAlloc(image->numCols,image->numRows,PS_TYPE_F32); - psMetadataAddPtr(data->OTA_solutions,PS_LIST_TAIL, chipName, PS_DATA_UNKNOWN | PS_META_REPLACE, "OTA solution element", solution); + psMetadataAddPtr(data->OTA_solutions,PS_LIST_TAIL, + chipName, PS_DATA_UNKNOWN | PS_META_REPLACE, + "OTA solution element", solution); } else { // We are not fitting OTAs, so read the one that should be saved on OTApath. @@ -203,6 +226,11 @@ psStringAppend(&solutionFileName, ".%s.fits",chipName); psFits *solutionFits = psFitsOpen(solutionFileName,"r"); - psImage *solution = psFitsReadImage(solutionFits,region,0); - psMetadataAddPtr(data->OTA_solutions,PS_LIST_TAIL, chipName, PS_DATA_UNKNOWN | PS_META_REPLACE, "OTA solution element", solution); + psImage *solution = psFitsReadImage(solutionFits,psRegionSet(0,0,0,0),0); + psMetadataAddImage(data->OTA_solutions,PS_LIST_TAIL, + chipName, 0, + "OTA solution element", solution); + psFree(solutionFits); + psFree(solutionFileName); + psFree(solution); } } @@ -210,17 +238,18 @@ // FIX Hopefully this isn't going to destroy everything. // Chip - if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { - psError(psErrorCodeLast(), false, "Error saving data to files."); - return false; - } +/* if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { */ +/* psError(psErrorCodeLast(), false, "Error saving data to files."); */ +/* return false; */ +/* } */ } // FPA - if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { - psError(psErrorCodeLast(), false, "Error saving data to files."); - return false; - } +/* if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { */ +/* psError(psErrorCodeLast(), false, "Error saving data to files."); */ +/* return false; */ +/* } */ // Add this smf to the smf array. psArrayAdd(data->smfs,data->smfs->n,psMemIncrRefCounter(smfFile)); + psFree(smfFile); } Index: /trunk/ppBackground/src/ppBackgroundStackData.c =================================================================== --- /trunk/ppBackground/src/ppBackgroundStackData.c (revision 36581) +++ /trunk/ppBackground/src/ppBackgroundStackData.c (revision 36582) @@ -21,4 +21,6 @@ psFree(data->modelMap); + + psFree(data->stack_data); psFree(data->stacks); psFree(data->outRoot); @@ -45,5 +47,6 @@ data->dec_min = 1e9; data->dec_max = -1e9; - + + data->stack_data = psArrayAlloc(0); data->stacks = psArrayAlloc(0); data->outRoot = NULL; Index: /trunk/ppBackground/src/ppBackgroundStackLoop.c =================================================================== --- /trunk/ppBackground/src/ppBackgroundStackLoop.c (revision 36581) +++ /trunk/ppBackground/src/ppBackgroundStackLoop.c (revision 36582) @@ -13,5 +13,5 @@ ) { - // pmConfig *config = data->config; // Configuration data + pmConfig *config = data->config; // Configuration data int i; // @@ -63,10 +63,114 @@ // Loop over the input images, and apply the models to construct the restored versions. - for (i = 0; i < data->stacks->n; i++) { + for (i = 0; i < data->stack_data->n; i++) { + pmFPAfile *stack = data->stack_data->data[i]; + pmFPAview *view = pmFPAviewAlloc(0); // Define output image - // Calculate model for each pixel of output - // Close output image + pmFPAfile *stack_model = pmFPAfileDefineOutput(config,stack->fpa,"PPBACKGROUND.STACK.MODEL"); + + if (!stack_model) { + psError(psErrorCodeLast(), false, "Unable to generate output model"); + return (false); + } + + pmFPAfile *stack_corr = pmFPAfileDefineOutput(config,stack->fpa,"PPBACKGROUND.STACK.OUTPUT"); + if (!stack_corr) { + psError(psErrorCodeLast(), false, "Unable to generate output result"); + return (false); + } + stack_model->save = true; + stack_corr->save = true; + + // Iterate over the images. + pmChip *chip; + while ((chip = pmFPAviewNextChip(view, stack->fpa, 1))) { + if (!chip->process || !chip->file_exists) { + continue; + } + + if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { + psError(psErrorCodeLast(), false, "load failure for Chip"); + return(false); + } + + pmCell *cell; + + while ((cell = pmFPAviewNextCell(view, stack->fpa, 1)) != NULL) { + psLogMsg ("ppImageLoop", 5, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process); + if (!cell->process || !cell->file_exists) { + continue; + } + if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { + psError(psErrorCodeLast(), false, "load failure for Cell"); + return(false); + } + if (cell->readouts->n > 1) { + psWarning ("Skipping Video Cell for ppImageDetrendReadout"); + continue; + } + + // process each of the readouts + pmReadout *readout; // Readout from cell + while ((readout = pmFPAviewNextReadout (view, stack->fpa, 1)) != NULL) { + if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { + psError(psErrorCodeLast(), false, "load failure for Readout"); + return(false); + } + if (!readout->data_exists) { + continue; + } + + pmReadout *model = pmFPAfileThisReadout(config->files,view,"PPBACKGROUND.STACK.MODEL"); + pmReadout *corr = pmFPAfileThisReadout(config->files,view,"PPBACKGROUND.STACK.OUTPUT"); + + psPlane *pix = psPlaneAlloc(); // Pixel coordinates on chip + psPlane *fp = psPlaneAlloc(); // Focal plane coordinates + psPlane *tp = psPlaneAlloc(); // Tangent plane coordinates + + int x,y; + for (y = 0; y < readout->image->numRows; y++) { + for (x = 0; x < readout->image->numCols; x++) { + // Calculate model for each pixel of output + psPlaneTransformApply(fp, chip->toFPA, pix); + psPlaneTransformApply(tp, stack->fpa->toTPA, fp); + + model->image->data.F32[y][x] = psImageMapEval(data->modelMap,tp->x,tp->y); + corr->image->data.F32[y][x] = readout->image->data.F32[y][x] + model->image->data.F32[y][x]; + } + } + psFree(pix); + psFree(fp); + psFree(tp); + } // Close readout + + // Close output image + // Close cells (XXX shouldn't pmFPAfileClose iterate down as needed?) + view->cell = -1; + while ((cell = pmFPAviewNextCell(view, stack->fpa, 1)) != NULL) { + if (!cell->process || !cell->file_exists) { + continue; + } + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + psError(psErrorCodeLast(), false, "save failure for Cell"); + return(false); + } + } + } // Close Cell + + // Close chip + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + psError(psErrorCodeLast(), false, "save failure for Chip"); + return(false); + } + } // Close chip. + // Output and Close FPA + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + psError(psErrorCodeLast(), false, "save failure for FPA"); + return(false); + } + psFree(view); } - + + return(true); }