Index: /branches/czw_branch/20100817/ippconfig/gpc1/ppImage.config =================================================================== --- /branches/czw_branch/20100817/ippconfig/gpc1/ppImage.config (revision 29485) +++ /branches/czw_branch/20100817/ippconfig/gpc1/ppImage.config (revision 29486) @@ -262,4 +262,6 @@ NOISEMAP METADATA END + LINEARITY METADATA + END END Index: /branches/czw_branch/20100817/ippconfig/recipes/filerules-mef.mdc =================================================================== --- /branches/czw_branch/20100817/ippconfig/recipes/filerules-mef.mdc (revision 29485) +++ /branches/czw_branch/20100817/ippconfig/recipes/filerules-mef.mdc (revision 29486) @@ -63,4 +63,5 @@ PPIMAGE.FRINGE INPUT @DETDB CHIP FRINGE PPIMAGE.SHUTTER INPUT @DETDB CHIP IMAGE +PPIMAGE.LINEARITY INPUT @DETDB CHIP LINEARITY ## Files used by ppMerge Index: /branches/czw_branch/20100817/ippconfig/recipes/filerules-simple.mdc =================================================================== --- /branches/czw_branch/20100817/ippconfig/recipes/filerules-simple.mdc (revision 29485) +++ /branches/czw_branch/20100817/ippconfig/recipes/filerules-simple.mdc (revision 29486) @@ -25,4 +25,5 @@ PPIMAGE.FRINGE INPUT @DETDB CHIP FRINGE PPIMAGE.SHUTTER INPUT @DETDB CHIP IMAGE +PPIMAGE.LINEARITY INPUT @DETDB CHIP LINEARITY ## Files used by ppMerge Index: /branches/czw_branch/20100817/ippconfig/recipes/filerules-split.mdc =================================================================== --- /branches/czw_branch/20100817/ippconfig/recipes/filerules-split.mdc (revision 29485) +++ /branches/czw_branch/20100817/ippconfig/recipes/filerules-split.mdc (revision 29486) @@ -36,4 +36,5 @@ PPIMAGE.FRINGE INPUT @DETDB CHIP FRINGE PPIMAGE.SHUTTER INPUT @DETDB CHIP IMAGE +PPIMAGE.LINEARITY INPUT @DETDB CHIP LINEARITY ## Files used by ppMerge Index: /branches/czw_branch/20100817/ppImage/src/ppImage.h =================================================================== --- /branches/czw_branch/20100817/ppImage/src/ppImage.h (revision 29485) +++ /branches/czw_branch/20100817/ppImage/src/ppImage.h (revision 29486) @@ -155,5 +155,6 @@ bool ppImageDetrendBias(pmReadout *inputReadout, pmReadout *bias, pmReadout *dark, ppImageOptions *options); -bool ppImageDetrendNonLinear(pmReadout *input, ppImageOptions *options); +//bool ppImageDetrendNonLinear(pmReadout *input, ppImageOptions *options); +bool ppImageDetrendNonLinear(pmReadout *input, pmFPAview *linearity, pmConfig *config); bool ppImageDetrendNonLinearLookup(pmReadout *input, psMetadataItem *dataItem); bool ppImageDetrendNonLinearPolynomial(pmReadout *input, psMetadataItem *dataItem); Index: /branches/czw_branch/20100817/ppImage/src/ppImageDefineFile.c =================================================================== --- /branches/czw_branch/20100817/ppImage/src/ppImageDefineFile.c (revision 29485) +++ /branches/czw_branch/20100817/ppImage/src/ppImageDefineFile.c (revision 29486) @@ -15,5 +15,5 @@ file = pmFPAfileDefineFromArgs(&status, config, filerule, argname); if (!status) { - psError(PS_ERR_UNKNOWN, false, "failed to load file definition"); + psError(PS_ERR_UNKNOWN, false, "failed to load file definition ARG LIST"); return false; } @@ -23,5 +23,5 @@ file = pmFPAfileDefineFromRun(&status, NULL, config, filerule); if (!status) { - psError(PS_ERR_UNKNOWN, false, "failed to load file definition"); + psError(PS_ERR_UNKNOWN, false, "failed to load file definition RUN"); return false; } @@ -31,5 +31,5 @@ file = pmFPAfileDefineFromConf(&status, config, filerule); if (!status) { - psError(PS_ERR_UNKNOWN, false, "failed to load file definition"); + psError(PS_ERR_UNKNOWN, false, "failed to load file definition CONFIG"); return false; } @@ -39,5 +39,5 @@ file = pmFPAfileDefineFromDetDB(&status, config, filerule, input, detrendType); if (!status) { - psError(PS_ERR_UNKNOWN, false, "failed to load file definition"); + psError(PS_ERR_UNKNOWN, false, "failed to load file definition DETREND"); return false; } Index: /branches/czw_branch/20100817/ppImage/src/ppImageDetrendFree.c =================================================================== --- /branches/czw_branch/20100817/ppImage/src/ppImageDetrendFree.c (revision 29485) +++ /branches/czw_branch/20100817/ppImage/src/ppImageDetrendFree.c (revision 29486) @@ -13,4 +13,5 @@ "PPIMAGE.FLAT", "PPIMAGE.SHUTTER", + "PPIMAGE.LINEARITY", NULL }; @@ -23,4 +24,5 @@ pmFPAfile *file = psMetadataLookupPtr(&status, config->files, detrendTypes[i]); // File of interest + psTrace("pmFPAfileFree",1,"Working on %s\n",detrendTypes[i]); if (!file) continue; // not all detrends are used in any given run Index: /branches/czw_branch/20100817/ppImage/src/ppImageDetrendNonLinear.c =================================================================== --- /branches/czw_branch/20100817/ppImage/src/ppImageDetrendNonLinear.c (revision 29485) +++ /branches/czw_branch/20100817/ppImage/src/ppImageDetrendNonLinear.c (revision 29486) @@ -4,4 +4,5 @@ #include "ppImage.h" + bool ppImageDetrendNonLinearPolynomial(pmReadout *input, psMetadataItem *dataItem) { @@ -44,71 +45,99 @@ } -bool ppImageDetrendNonLinear(pmReadout *input, ppImageOptions *options) { - psMetadataItem *concept; - pmCell *cell = input->parent; +bool ppImageDetrendNonLinear(pmReadout *input, pmFPAview *detview, pmConfig *config) { + bool status; - switch (options->nonLinearType) { - case PS_DATA_VECTOR: - ppImageDetrendNonLinearPolynomial (input, options->nonLinearData); - return true; + pmFPAfile *linearity_file = psMetadataLookupPtr(&status,config->files,"PPIMAGE.LINEARITY"); + psFits *linearity_fits = linearity_file->fits; - case PS_DATA_STRING: - ppImageDetrendNonLinearLookup (input, options->nonLinearData); - return true; + if (!psFitsMoveExtName(linearity_fits,psMetadataLookupStr(&status,input->parent->concepts,"CELL.NAME"))) { + psError(PS_ERR_IO, false, "Unable to move to non-linearity table %s",psMetadataLookupStr(&status,input->parent->concepts,"CELL.NAME")); + return(false); + } + + psArray *table = psFitsReadTable(linearity_fits); + if (!table) { + psError(PS_ERR_IO, false, "Unable to read non-linearity table.\n"); + return(false); + } - case PS_DATA_METADATA: - // XXX EAM: this is somewhat confusing : let's wrap in a function when i understand it + // It might be better to pack lookup table here... + // Why? I only use that lookup table once for the single cell it matches. + + if (!pmNonLinearityApply(input,table)) { + psError(PS_ERR_UNKNOWN, false, "Unable to apply non-linearity corrections.\n"); + return(false); + } - // Go looking for the value in the hierarchy - concept = psMetadataLookup(cell->concepts, options->nonLinearSource); - if (! concept) { - pmChip *chip = cell->parent;// Parent chip - concept = psMetadataLookup(chip->concepts, options->nonLinearSource); - if (! concept) { - pmFPA *fpa = chip->parent; // Parent FPA - concept = psMetadataLookup(fpa->concepts, options->nonLinearSource); - if (! concept) { - psLogMsg("phase2", PS_LOG_WARN, "Unable to find value of concept %s " - "for non-linearity correction --- ignored.\n", (char *)options->nonLinearSource); - return false; - } - } - } - if (concept->type != PS_DATA_STRING) { - psLogMsg("phase2", PS_LOG_WARN, "Type for concept %p isn't STRING, as" - " expected for non-linearity correction --- ignored.\n", - concept); - return false; - } + return(true); +} + /* bool ppImageDetrendNonLinear(pmReadout *input, ppImageOptions *options) { */ +/* psMetadataItem *concept; */ +/* pmCell *cell = input->parent; */ - // Get the value of the concept - psString conceptValue = concept->data.V; - psMetadata *folder = (psMetadata *)options->nonLinearData->data.V; - psMetadataItem *optionItem = psMetadataLookup(folder, conceptValue); - if (!optionItem) { - psLogMsg("phase2", PS_LOG_WARN, "Unable to find %s in NONLIN.DATA" - " --- ignored.\n", conceptValue); - return false; - } +/* switch (options->nonLinearType) { */ +/* case PS_DATA_VECTOR: */ +/* ppImageDetrendNonLinearPolynomial (input, options->nonLinearData); */ +/* return true; */ - switch (optionItem->type) { - case PS_DATA_VECTOR: - ppImageDetrendNonLinearPolynomial (input, optionItem); - return true; - case PS_DATA_STRING: - ppImageDetrendNonLinearLookup (input, optionItem); - return true; - default: - psLogMsg("phase2", PS_LOG_WARN, "Non-linearity correction " - "desired but unable to interpret NONLIN.DATA for %s" - " --- ignored\n", conceptValue); - return false; - } - default: - psAbort("Invalid options->nonLinearType"); - } - return true; -} +/* case PS_DATA_STRING: */ +/* ppImageDetrendNonLinearLookup (input, options->nonLinearData); */ +/* return true; */ +/* case PS_DATA_METADATA: */ +/* // XXX EAM: this is somewhat confusing : let's wrap in a function when i understand it */ + +/* // Go looking for the value in the hierarchy */ +/* concept = psMetadataLookup(cell->concepts, options->nonLinearSource); */ +/* if (! concept) { */ +/* pmChip *chip = cell->parent;// Parent chip */ +/* concept = psMetadataLookup(chip->concepts, options->nonLinearSource); */ +/* if (! concept) { */ +/* pmFPA *fpa = chip->parent; // Parent FPA */ +/* concept = psMetadataLookup(fpa->concepts, options->nonLinearSource); */ +/* if (! concept) { */ +/* psLogMsg("phase2", PS_LOG_WARN, "Unable to find value of concept %s " */ +/* "for non-linearity correction --- ignored.\n", (char *)options->nonLinearSource); */ +/* return false; */ +/* } */ +/* } */ +/* } */ + +/* if (concept->type != PS_DATA_STRING) { */ +/* psLogMsg("phase2", PS_LOG_WARN, "Type for concept %p isn't STRING, as" */ +/* " expected for non-linearity correction --- ignored.\n", */ +/* concept); */ +/* return false; */ +/* } */ + +/* // Get the value of the concept */ +/* psString conceptValue = concept->data.V; */ +/* psMetadata *folder = (psMetadata *)options->nonLinearData->data.V; */ +/* psMetadataItem *optionItem = psMetadataLookup(folder, conceptValue); */ +/* if (!optionItem) { */ +/* psLogMsg("phase2", PS_LOG_WARN, "Unable to find %s in NONLIN.DATA" */ +/* " --- ignored.\n", conceptValue); */ +/* return false; */ +/* } */ + +/* switch (optionItem->type) { */ +/* case PS_DATA_VECTOR: */ +/* ppImageDetrendNonLinearPolynomial (input, optionItem); */ +/* return true; */ +/* case PS_DATA_STRING: */ +/* ppImageDetrendNonLinearLookup (input, optionItem); */ +/* return true; */ +/* default: */ +/* psLogMsg("phase2", PS_LOG_WARN, "Non-linearity correction " */ +/* "desired but unable to interpret NONLIN.DATA for %s" */ +/* " --- ignored\n", conceptValue); */ +/* return false; */ +/* } */ +/* default: */ +/* psAbort("Invalid options->nonLinearType"); */ +/* } */ +/* return true; */ +/* } */ + Index: /branches/czw_branch/20100817/ppImage/src/ppImageDetrendReadout.c =================================================================== --- /branches/czw_branch/20100817/ppImage/src/ppImageDetrendReadout.c (revision 29485) +++ /branches/czw_branch/20100817/ppImage/src/ppImageDetrendReadout.c (revision 29486) @@ -51,11 +51,14 @@ } - -# if 0 // Non-linearity correction if (options->doNonLin) { - ppImageDetrendNonLinear(detrend->input, input, options); - } -# endif + // linearity = pmFPAfileThisReadout(config->files, detview, "PPIMAGE.LINEARITY"); + if (!ppImageDetrendNonLinear(input,detview,config)) { + psError(PS_ERR_UNKNOWN, false, "Unable to correct NonLinearity"); + psFree(detview); + return(false); + } +/* ppImageDetrendNonLinear(detrend->input, input, options); */ + } // set up the dark and bias Index: /branches/czw_branch/20100817/ppImage/src/ppImageLoop.c =================================================================== --- /branches/czw_branch/20100817/ppImage/src/ppImageLoop.c (revision 29485) +++ /branches/czw_branch/20100817/ppImage/src/ppImageLoop.c (revision 29486) @@ -152,5 +152,5 @@ ESCAPE("Unable to free detrend images"); } - + // Apply the fringe correction if (options->doFringe) { Index: /branches/czw_branch/20100817/ppImage/src/ppImageParseCamera.c =================================================================== --- /branches/czw_branch/20100817/ppImage/src/ppImageParseCamera.c (revision 29485) +++ /branches/czw_branch/20100817/ppImage/src/ppImageParseCamera.c (revision 29486) @@ -73,4 +73,13 @@ return NULL; } + } + + if (options->doNonLin) { + if (!ppImageDefineFile(config, input->fpa, "PPIMAGE.LINEARITY", "LINEARITY", + PM_FPA_FILE_LINEARITY, PM_DETREND_TYPE_LINEARITY)) { + psError(PS_ERR_IO, false, "Can't find a non-linearity correction source"); + psFree(options); + return NULL; + } } Index: /branches/czw_branch/20100817/psLib/src/types/psMetadata.c =================================================================== --- /branches/czw_branch/20100817/psLib/src/types/psMetadata.c (revision 29485) +++ /branches/czw_branch/20100817/psLib/src/types/psMetadata.c (revision 29486) @@ -109,5 +109,5 @@ type = metadataItem->type; - + // fprintf(stderr,"%s %s %d\n",metadataItem->name,metadataItem->comment,type); psMemDecrRefCounter(metadataItem->name); psMemDecrRefCounter(metadataItem->comment); Index: /branches/czw_branch/20100817/psModules/src/camera/pmFPAfile.c =================================================================== --- /branches/czw_branch/20100817/psModules/src/camera/pmFPAfile.c (revision 29485) +++ /branches/czw_branch/20100817/psModules/src/camera/pmFPAfile.c (revision 29486) @@ -36,5 +36,5 @@ return; } - + psTrace ("pmFPAfileFree", 5, "freeing %s %ld\n", file->name,(psU64) file->fits); psAssert(!fpaFileFreeStrict || file->fits == NULL, "File %s wasn't closed.", file->name); @@ -523,4 +523,7 @@ if (!strcasecmp(type, "HEADER")) { return PM_FPA_FILE_HEADER; + } + if (!strcasecmp(type, "LINEARITY")) { + return PM_FPA_FILE_LINEARITY; } // deprecate this? Index: /branches/czw_branch/20100817/psModules/src/camera/pmFPAfile.h =================================================================== --- /branches/czw_branch/20100817/psModules/src/camera/pmFPAfile.h (revision 29485) +++ /branches/czw_branch/20100817/psModules/src/camera/pmFPAfile.h (revision 29486) @@ -50,4 +50,5 @@ PM_FPA_FILE_SRCTEXT, PM_FPA_FILE_PATTERN, + PM_FPA_FILE_LINEARITY, } pmFPAfileType; Index: /branches/czw_branch/20100817/psModules/src/camera/pmFPAfileFitsIO.c =================================================================== --- /branches/czw_branch/20100817/psModules/src/camera/pmFPAfileFitsIO.c (revision 29485) +++ /branches/czw_branch/20100817/psModules/src/camera/pmFPAfileFitsIO.c (revision 29486) @@ -145,4 +145,5 @@ case PM_FPA_FILE_FRINGE: case PM_FPA_FILE_DARK: + case PM_FPA_FILE_LINEARITY: case PM_FPA_FILE_CMP: case PM_FPA_FILE_CMF: Index: /branches/czw_branch/20100817/psModules/src/camera/pmFPAfileIO.c =================================================================== --- /branches/czw_branch/20100817/psModules/src/camera/pmFPAfileIO.c (revision 29485) +++ /branches/czw_branch/20100817/psModules/src/camera/pmFPAfileIO.c (revision 29486) @@ -224,4 +224,5 @@ case PM_FPA_FILE_JPEG: case PM_FPA_FILE_KAPA: + case PM_FPA_FILE_LINEARITY: break; default: @@ -551,4 +552,5 @@ case PM_FPA_FILE_ASTROM_MODEL: case PM_FPA_FILE_ASTROM_REFSTARS: + case PM_FPA_FILE_LINEARITY: psTrace ("psModules.camera", 5, "closing %s (%s) (%d:%d:%d)\n", file->filename, file->name, view->chip, view->cell, view->readout); status = psFitsClose (file->fits); @@ -566,4 +568,5 @@ case PM_FPA_FILE_JPEG: case PM_FPA_FILE_KAPA: + break; @@ -610,4 +613,5 @@ case PM_FPA_FILE_FRINGE: case PM_FPA_FILE_DARK: + // status = pmFPAviewFreeData(view, file); break; @@ -627,4 +631,5 @@ case PM_FPA_FILE_JPEG: case PM_FPA_FILE_KAPA: + case PM_FPA_FILE_LINEARITY: psTrace ("psModules.camera", 5, "nothing to free for %s (%s)\n", file->filename, file->name); return true; @@ -782,4 +787,5 @@ case PM_FPA_FILE_ASTROM_MODEL: case PM_FPA_FILE_ASTROM_REFSTARS: + case PM_FPA_FILE_LINEARITY: psTrace ("psModules.camera", 5, "opening %s (%s) (%d:%d:%d)\n", file->filename, file->name, view->chip, view->cell, view->readout); Index: /branches/czw_branch/20100817/psModules/src/detrend/pmDetrendDB.c =================================================================== --- /branches/czw_branch/20100817/psModules/src/detrend/pmDetrendDB.c (revision 29485) +++ /branches/czw_branch/20100817/psModules/src/detrend/pmDetrendDB.c (revision 29486) @@ -106,4 +106,5 @@ DETREND_STRING_CASE(ASTROM); DETREND_STRING_CASE(NOISEMAP); + DETREND_STRING_CASE(LINEARITY); default: return NULL; Index: /branches/czw_branch/20100817/psModules/src/detrend/pmDetrendDB.h =================================================================== --- /branches/czw_branch/20100817/psModules/src/detrend/pmDetrendDB.h (revision 29485) +++ /branches/czw_branch/20100817/psModules/src/detrend/pmDetrendDB.h (revision 29486) @@ -36,4 +36,5 @@ PM_DETREND_TYPE_ASTROM, PM_DETREND_TYPE_NOISEMAP, + PM_DETREND_TYPE_LINEARITY, } pmDetrendType; Index: /branches/czw_branch/20100817/psModules/src/detrend/pmNonLinear.c =================================================================== --- /branches/czw_branch/20100817/psModules/src/detrend/pmNonLinear.c (revision 29485) +++ /branches/czw_branch/20100817/psModules/src/detrend/pmNonLinear.c (revision 29486) @@ -45,4 +45,28 @@ } } + +# define PS_BIN_INTERPOLATE(RESULT, VECTOR, BOUNDS, BIN, VALUE) { \ + float dX, dY, Xo, Yo, Xt; \ + if (BIN == BOUNDS->n - 1) { \ + dX = 0.5*(BOUNDS->data.F32[BIN+1] - BOUNDS->data.F32[BIN-1]); \ + dY = VECTOR->data.F32[BIN] - VECTOR->data.F32[BIN-1]; \ + Xo = 0.5*(BOUNDS->data.F32[BIN+1] + BOUNDS->data.F32[BIN]); \ + Yo = VECTOR->data.F32[BIN]; \ + } else { \ + dX = 0.5*(BOUNDS->data.F32[BIN+2] - BOUNDS->data.F32[BIN]); \ + dY = VECTOR->data.F32[BIN+1] - VECTOR->data.F32[BIN]; \ + Xo = 0.5*(BOUNDS->data.F32[BIN+1] + BOUNDS->data.F32[BIN]); \ + Yo = VECTOR->data.F32[BIN]; \ + } \ + if (dY != 0) { \ + Xt = (VALUE - Yo)*dX/dY + Xo; \ + } else { \ + Xt = Xo; \ + } \ + Xt = PS_MIN (BOUNDS->data.F32[BIN+1], PS_MAX(BOUNDS->data.F32[BIN], Xt)); \ + psTrace("pmNonLinear", 6, "(Xo, Yo, dX, dY, Xt, Yt) is (%.2f %.2f %.2f %.2f %.2f %.2f)\n", \ + Xo, Yo, dX, dY, Xt, VALUE); \ + RESULT = Xt; } + pmReadout *pmNonLinearityLookup(pmReadout *inputReadout, const psVector *inFlux, const psVector *outFlux) { @@ -92,2 +116,182 @@ return inputReadout; } + +bool pmNonLinearityApply(pmReadout *inputReadout, psArray *Ltab) +{ + PS_ASSERT_PTR_NON_NULL(inputReadout, false); + PS_ASSERT_PTR_NON_NULL(inputReadout->image, false); + PS_ASSERT_IMAGE_TYPE(inputReadout->image, PS_TYPE_F32, false); + PS_ASSERT_PTR_NON_NULL(Ltab, false); + + // psS32 tableSize = Ltab->n; + + psImage *image = inputReadout->image; + psS32 numSamples = 39; + + psS32 numBorder = 10; + // Load default data. + psVector *default_row_correction_fluxes = psVectorAlloc(numSamples,PS_TYPE_F32); + psVector *default_col_correction_fluxes = psVectorAlloc(numSamples,PS_TYPE_F32); + psVector *row_correction_fluxes; + psVector *col_correction_fluxes; + + psVector *default_row_correction_factors = psVectorAlloc(numSamples,PS_TYPE_F32); + psVector *default_col_correction_factors = psVectorAlloc(numSamples,PS_TYPE_F32); + psVector *row_correction_factors; + psVector *col_correction_factors; + + int n = 0; + int m = 0; + for (int k = 0; k < Ltab->n; k++) { // Begin load default tables + psMetadata *row = Ltab->data[k]; + if (psMetadataLookupS32(NULL,row,"POSITION") != -1) { + continue; + } + if (psMetadataLookupS32(NULL,row,"DIRECTION") == 0) { + psVectorSet(default_row_correction_fluxes,n,psMetadataLookupF32(NULL,row,"FLUX")); + psVectorSet(default_row_correction_factors,n,psMetadataLookupF32(NULL,row,"FACTOR")); + n++; + } + else { + psVectorSet(default_col_correction_fluxes,m,psMetadataLookupF32(NULL,row,"FLUX")); + psVectorSet(default_col_correction_factors,m,psMetadataLookupF32(NULL,row,"FACTOR")); + m++; + } + } // End load default tables + + for (int i = 0; i < image->numRows; i++) { // Loop over rows : note: problem with discontinuity here + n = 0; + if ((i < numBorder)||(image->numRows - i < numBorder)) { // Load row correction data + row_correction_fluxes = psVectorAlloc(numSamples,PS_TYPE_F32); + row_correction_factors = psVectorAlloc(numSamples,PS_TYPE_F32); + + for (int k = 0; k < Ltab->n; k++) { + psMetadata *row = Ltab->data[k]; + if ((psMetadataLookupS32(NULL,row,"DIRECTION") != 0)|| + (psMetadataLookupS32(NULL,row,"POSITION") != i + 1)) { // image data is 0 indexed, but correction data is 1 indexed. + continue; + } + psVectorSet(row_correction_fluxes,n,psMetadataLookupF32(NULL,row,"FLUX")); + psVectorSet(row_correction_factors,n,psMetadataLookupF32(NULL,row,"FACTOR")); + // psTrace("psModules.camera",6,"rows: %d\n",n); + n++; + } + } // End load specific row + else { + row_correction_fluxes = default_row_correction_fluxes; + row_correction_factors = default_row_correction_factors; + } // End load default row + + for (int j = 0; j < image->numCols; j++) { // Loop over columns + m = 0; + + if ((j < numBorder)||(image->numCols - j < numBorder)) { // Load col correction data + col_correction_fluxes = psVectorAlloc(numSamples,PS_TYPE_F32); + col_correction_factors = psVectorAlloc(numSamples,PS_TYPE_F32); + for (int k = 0; k < Ltab->n; k++) { + psMetadata *row = Ltab->data[k]; + if ((psMetadataLookupS32(NULL,row,"DIRECTION") != 1)|| + (psMetadataLookupS32(NULL,row,"POSITION") != j + 1)) { + continue; + } + psVectorSet(col_correction_fluxes,m,psMetadataLookupF32(NULL,row,"FLUX")); + psVectorSet(col_correction_factors,m,psMetadataLookupF32(NULL,row,"FACTOR")); + // psTrace("psModules.camera",6,"columns: %d\n",m); + m++; + } + } // End load specific col + else { + col_correction_fluxes = default_col_correction_fluxes; + col_correction_factors = default_col_correction_factors; + } // End load default col + + // Calculate correction factor contribution for this pixel. + psF32 factor_row = pmNonLinearityMeasure(image->data.F32[i][j], + row_correction_fluxes,row_correction_factors); + psF32 factor_col = pmNonLinearityMeasure(image->data.F32[i][j], + col_correction_fluxes,col_correction_factors); + + if (((i <= 9)&&(j <= 9))||((i == 9)&&(j == 5))) { // Print out if we're looking at a test case. + psTrace("psModules.camera",6,"Linearity: %d %d %s %f %f %f %d %d\n",i,j, + psMetadataLookupStr(NULL,inputReadout->parent->concepts,"CELL.NAME"), + image->data.F32[i][j],factor_row,factor_col,numBorder,numSamples); + psTrace("psModules.camera",6,"Linearity: R: %d %d %d C: %d %d %d\n", + i,(i < numBorder),(image->numRows - i), + j,(j < numBorder),(image->numCols - j)); + +/* psTrace("psModules.camera",6,"Linearity: V: "); */ +/* for (int k = 0; k < numSamples; k++) { */ +/* psTrace("psModules.camera",6,"(%f %f) ", */ +/* col_correction_fluxes->data.F32[k],col_correction_factors->data.F32[k]); */ +/* } */ +/* psTrace("psModules.camera",6,"\n"); */ + } // End Test case + + // Apply correction to image data + image->data.F32[i][j] = image->data.F32[i][j] * ( factor_row * factor_col ) ; + + // Free correction if we allocated + if ((j < numBorder)||(image->numCols - j < numBorder)) { + psFree(col_correction_fluxes); + psFree(col_correction_factors); + } + + } // End loop over columns + + // Free correction if we allocated + if ((i < numBorder)||(image->numRows - i < numBorder)) { + psFree(row_correction_fluxes); + psFree(row_correction_factors); + } + + } // End loop over rows + + psFree(default_row_correction_fluxes); + psFree(default_row_correction_factors); + psFree(default_col_correction_fluxes); + psFree(default_col_correction_factors); + + return(true); +} + +psF32 pmNonLinearityMeasure(psF32 flux, psVector *correction_fluxes, psVector *correction_factors) { + // psS32 numPixels = 0; + psF32 result = 0; + psU32 bin = 0; + + bin = correction_fluxes->n - 1; + if (flux < correction_fluxes->data.F32[0]) { + return(1.0); + } + if (flux > correction_fluxes->data.F32[bin]) { + return(1.0); + } + + for (int i = 0; i < correction_fluxes->n - 1; i++) { + if ((flux >= correction_fluxes->data.F32[i])&& + (flux < correction_fluxes->data.F32[i+1])) { + result = correction_factors->data.F32[i] + + (flux - correction_fluxes->data.F32[i]) * + ((correction_factors->data.F32[i+1] - correction_factors->data.F32[i]) / + (correction_fluxes->data.F32[i+1] - correction_fluxes->data.F32[i])); + continue; + } + } + +/* PS_BIN_FOR_VALUE(bin,correction_fluxes,flux); */ +/* if ((bin < 0)||(bin > correction_fluxes->n)) { */ +/* return(1.0); */ +/* } */ + +/* PS_BIN_INTERPOLATE(result,correction_fluxes,correction_factors,bin,flux); */ + if (!isfinite(result)) { + result = 1.0; + } + if (result <= 0) { + result = 1.0; + } + return(result); +} + + + Index: /branches/czw_branch/20100817/psModules/src/detrend/pmNonLinear.h =================================================================== --- /branches/czw_branch/20100817/psModules/src/detrend/pmNonLinear.h (revision 29485) +++ /branches/czw_branch/20100817/psModules/src/detrend/pmNonLinear.h (revision 29486) @@ -31,5 +31,6 @@ const psVector *outFlux ///< Table column with output fluxes ); - +bool pmNonLinearityApply(pmReadout *inputReadout, psArray *Ltab); +psF32 pmNonLinearityMeasure(psF32 flux, psVector *correction_fluxes, psVector *correction_factors); /// @} #endif