Index: trunk/psastro/src/psastroMaskUpdates.c =================================================================== --- trunk/psastro/src/psastroMaskUpdates.c (revision 26897) +++ trunk/psastro/src/psastroMaskUpdates.c (revision 28043) @@ -19,8 +19,51 @@ } + +int colStart(psU16 corner_list,int i,int x_0,int y_0,int R,int max) { + if ((corner_list == 0x04)||(corner_list == 0x05)) { + return(0); + } + else { + return((int) fabs(x_0 - sqrt(pow(R,2) - pow(y_0 - i,2)))); + } +} +int colEnd(psU16 corner_list,int i,int x_0,int y_0,int R,int max) { + if ((corner_list == 0x04)||(corner_list == 0x05)) { + return((int) fabs(x_0 - sqrt(pow(R,2) - pow(y_0 - i,2)))); + } + else { + return(max); + } + +} +int rowStart(psU16 corner_list,int j,int x_0,int y_0,int R,int max) { + if ((corner_list == 0x01)||(corner_list == 0x02)|| + (corner_list == 0x03)||(corner_list == 0x0b)|| + (corner_list == 0x07)||(corner_list == 0x0f)) { + return(0); + } + else { + return((int) fabs(y_0 - sqrt(pow(R,2) - pow(x_0 - j,2)))); + } +} +int rowEnd(psU16 corner_list,int j,int x_0,int y_0,int R, int max) { + if (corner_list == 0x0f) { + return(max); + } + else if ((corner_list == 0x01)||(corner_list == 0x02)|| + (corner_list == 0x03)||(corner_list == 0x0b)|| + (corner_list == 0x07)) { + return((int) fabs(y_0 - sqrt(pow(R,2) - pow(x_0 - j,2)))); + } + else { + return(max); + } + +} + /** * create a mask or mask regions based on the collection of reference stars that * are in the vicinity of each chip */ -bool psastroMaskUpdates (pmConfig *config) { +bool psastroMaskUpdates (pmConfig *config, psMetadata *stats) { bool status; @@ -88,4 +131,18 @@ double REFSTAR_MASK_CROSSTALK_MAG_SLOPE= psMetadataLookupF32 (&status, recipe, "REFSTAR_MASK_CROSSTALK_MAG_SLOPE"); + // Mask stats variables + psU16 staticMaskVal = psMetadataLookupU32(&status, recipe, "MASKSTAT.STATIC"); + psU16 magicMaskVal = psMetadataLookupU32(&status, recipe, "MASKSTAT.MAGIC"); + psU16 dynamicMaskVal = psMetadataLookupU32(&status, recipe, "MASKSTAT.DYNAMIC"); + psU16 advisoryMaskVal = psMetadataLookupU32(&status, recipe, "MASKSTAT.ADVISORY"); + + psS32 Npix_valid = 0; + psS32 Npix_static = 0; + psS32 Npix_magic = 0; + psS32 Npix_dynamic = 0; + psS32 Npix_advisory = 0; + + psU16 corner_list = 0x00; + // select the input data sources pmFPAfile *input = psMetadataLookupPtr (NULL, config->files, "PSASTRO.INPUT"); @@ -119,4 +176,12 @@ pmFPA *fpaMask = outMask->fpa; + // Get camera specific mask stat options + psF32 FOV = psMetadataLookupF32(&status, fpaMask->camera, "FOV"); // Almost certainly doesn't exist. Placeholder for now. + psS32 CONSTANT_BLANK = psMetadataLookupS32(&status, fpaMask->camera, "NPIX_INTERCHIP"); // Need to add this as well. + + Npix_valid += CONSTANT_BLANK; + Npix_static += CONSTANT_BLANK; + + // select the reference mask fpa :: we use this to determine cell boundaries pmFPAfile *refMask = psMetadataLookupPtr (NULL, config->files, "PSASTRO.REFMASK"); @@ -363,6 +428,115 @@ } } - } - + } + + // Do the mask stats bit from cell down. + if (psMetadataLookupBool(&status,recipe,"MASK.STATS")) { + if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) ESCAPE; + while ((cell = pmFPAviewNextCell(view, fpa, 1)) != NULL) { + if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) ESCAPE; + psWarning("In CELL: FOV: %f CONSTANT: %d Mask bits: %d %d %d %d\n", + FOV,CONSTANT_BLANK,staticMaskVal,magicMaskVal,dynamicMaskVal,advisoryMaskVal); + + if (!cell->process || !cell->file_exists) {continue; } + while ((readout = pmFPAviewNextReadout(view, fpa, 1)) != NULL) { + pmReadout *readoutMask = pmFPAviewThisReadout (view, outMask->fpa); + if (!pmFPAfileIOChecks (config, viewMask, PM_FPA_BEFORE)) ESCAPE; + psWarning("In READOUT: FOV: %f CONSTANT: %d Mask bits: %d %d %d %d\n", + FOV,CONSTANT_BLANK,staticMaskVal,magicMaskVal,dynamicMaskVal,advisoryMaskVal); + + if (!readoutMask->data_exists) {continue; } + psPlane coordFPA; + psPlane coordCell; + + psImage *mask = readoutMask->mask; + if (!mask) {continue;} + psImageMaskType **maskData = mask->data.PS_TYPE_IMAGE_MASK_DATA; + // Dance coordinates around + // Calculate which corners fall within the field of view. + // 0x04 0x08 + // 0x01 0x02 + + coordFPA.x = 0.0; + coordFPA.y = 0.0; + psPlaneTransformApply(&coordCell,chip->fromFPA,&coordFPA); + if (pow(coordFPA.x,2) + pow(coordFPA.y,2) <= pow(FOV,2)) { + corner_list = corner_list | 0x01; + } + + coordCell.x = mask->numCols - 1; + psPlaneTransformApply(&coordCell,chip->fromFPA,&coordFPA); + if (pow(coordFPA.x,2) + pow(coordFPA.y,2) <= pow(FOV,2)) { + corner_list = corner_list | 0x02; + } + + coordCell.y = mask->numRows - 1; + psPlaneTransformApply(&coordCell,chip->fromFPA,&coordFPA); + if (pow(coordFPA.x,2) + pow(coordFPA.y,2) <= pow(FOV,2)) { + corner_list = corner_list | 0x08; + } + + coordCell.x = 0; + psPlaneTransformApply(&coordCell,chip->fromFPA,&coordFPA); + if (pow(coordFPA.x,2) + pow(coordFPA.y,2) <= pow(FOV,2)) { + corner_list = corner_list | 0x04; + } + + int x_0 = coordFPA.x; + int y_0 = coordFPA.y; + + + // Scan over the valid regions of the image and count masked pixels + if ((corner_list == 0x04)||(corner_list == 0x08)|| + (corner_list == 0x05)||(corner_list == 0x0a)) { + for (int i = 0; i < mask->numRows - 1; i++) { + for (int j = colStart(corner_list,i,x_0,y_0,FOV,mask->numCols - 1); + j < colEnd(corner_list,i,x_0,y_0,FOV,mask->numCols - 1); j++) { + Npix_valid++; + if (maskData[i][j] & staticMaskVal) { + Npix_static++; + continue; + } + if (maskData[i][j] & dynamicMaskVal) { + Npix_dynamic++; + continue; + } + if (maskData[i][j] & magicMaskVal) { + Npix_magic++; + continue; + } + if (maskData[i][j] & advisoryMaskVal) { + Npix_advisory++; + continue; + } + } + } + } + else { + for (int j = 0; j < mask->numCols - 1; j++) { + for (int i = rowStart(corner_list,j,x_0,y_0,FOV,mask->numRows - 1); + i < rowEnd(corner_list,j,x_0,y_0,FOV,mask->numRows - 1); i++) { + Npix_valid++; + if (maskData[i][j] & staticMaskVal) { + Npix_static++; + continue; + } + if (maskData[i][j] & dynamicMaskVal) { + Npix_dynamic++; + continue; + } + if (maskData[i][j] & magicMaskVal) { + Npix_magic++; + continue; + } + if (maskData[i][j] & advisoryMaskVal) { + Npix_advisory++; + continue; + } + } + } + } + } + } + } // output sequence for mask corresponding to this chip (XXX this may not be needed...) *viewMask = *view; @@ -392,4 +566,15 @@ } + psMetadataAddS32(stats,PS_LIST_TAIL, "MASKFRAC_NPIX", 0, + "Number of valid pixels", Npix_valid); + psMetadataAddF32(stats,PS_LIST_TAIL, "MASKFRAC_STATIC", 0, + "Fraction of pixels statically masked", (float) Npix_static / Npix_valid); + psMetadataAddF32(stats,PS_LIST_TAIL, "MASKFRAC_DYNAMIC", 0, + "Fraction of pixels dynamically masked", (float) Npix_dynamic / Npix_valid); + psMetadataAddF32(stats,PS_LIST_TAIL, "MASKFRAC_MAGIC", 0, + "Fraction of pixels magically masked", (float) Npix_magic / Npix_valid); + psMetadataAddF32(stats,PS_LIST_TAIL, "MASKFRAC_ADVISORY", 0, + "Fraction of pixels masked as an advisory", (float) Npix_advisory / Npix_valid); + // deactivate all files pmFPAfileActivate (config->files, false, NULL);