IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 29215


Ignore:
Timestamp:
Sep 23, 2010, 7:02:45 AM (16 years ago)
Author:
eugene
Message:

replace pmSubtractionMaskInvalid with pmReadoutMaskInvalid (takes test and set masks); change stack mask to propagate good (or suspect) and bad pixel masks; replace test blocks with CHECKPIX test macro; replace robust median for window with sample median (less prone to fail)

Location:
branches/eam_branches/ipp-20100823/psModules/src
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • branches/eam_branches/ipp-20100823/psModules/src/camera/pmFPAMaskWeight.c

    r26893 r29215  
    476476}
    477477
    478 
    479 
     478// find any pixels which are not already masked (with maskTest) which are not valid and raise maskSet bits
     479bool pmReadoutMaskInvalid (const pmReadout *readout, psImageMaskType maskTest, psImageMaskType maskSet) {
     480
     481    if (!readout) return true;
     482
     483    psImage *image = readout->image;
     484    psImage *mask  = readout->mask;
     485    psImage *variance = readout->variance;
     486    for (int y = 0; y < image->numRows; y++) {
     487        for (int x = 0; x < image->numCols; x++) {
     488            if (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskTest) continue;
     489            bool valid = false;
     490            valid = isfinite(image->data.F32[y][x]);
     491            if (valid && variance) {
     492                valid &= isfinite(variance->data.F32[y][x]);
     493            }
     494            if (valid) continue;
     495            mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] |= maskSet;
     496        }
     497    }
     498
     499    return true;
     500}
     501
     502// raise maskVal for any invalid pixels
    480503bool pmReadoutMaskNonfinite(pmReadout *readout, psImageMaskType maskVal)
    481504{
  • branches/eam_branches/ipp-20100823/psModules/src/camera/pmFPAMaskWeight.h

    r26076 r29215  
    9999    );
    100100
     101// find any pixels which are not already masked (with maskTest) which are not valid and raise maskSet bits
     102bool pmReadoutMaskInvalid (const pmReadout *readout, psImageMaskType maskTest, psImageMaskType maskSet);
     103
    101104/// Apply a mask to the image and variance map
    102105///
  • branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmStack.c

    r27427 r29215  
    1919
    2020#include <stdio.h>
     21#include <stdarg.h>
    2122#include <pslib.h>
    2223
     
    3435
    3536
    36 //#define TESTING                         // Enable test output
    37 //#define TEST_X 843-1                     // x coordinate to examine
    38 //#define TEST_Y 813-1                     // y coordinate to examine
    39 //#define TEST_RADIUS 0                    // Radius to examine
    40 
     37#define TESTING                         // Enable test output
     38// #define TEST_X 4963                       // x coordinate to examine
     39// #define TEST_Y 5353                       // y coordinate to examine
     40#define TEST_X 40                       // x coordinate to examine
     41#define TEST_Y 40                       // y coordinate to examine
     42#define TEST_RADIUS 0.5                 // Radius to examine
     43
     44# ifdef TESTING
     45# define CHECKPIX(XPIX,YPIX,MSG,...) { if (PS_SQR(XPIX - TEST_X) + PS_SQR(YPIX - TEST_Y) <= PS_SQR(TEST_RADIUS)) { fprintf(stderr,MSG,__VA_ARGS__); } }
     46# else
     47# define CHECKPIX(XPIX,YPIX,MSG,...) { }
     48# endif   
    4149
    4250// Data structure for use as a buffer when combining pixels
     
    250258                                      )
    251259{
    252 #ifdef TESTING
    253     if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    254         fprintf(stderr, "Marking image %d, pixel %d,%d for inspection\n", source, x, y);
    255     }
    256 #endif
     260    CHECKPIX(x, y, "Marking image %d, pixel %d,%d for inspection\n", source, x, y);
    257261    pmStackData *data = inputs->data[source]; // Stack data of interest
    258262    if (!data) {
     
    272276                                     )
    273277{
    274 #ifdef TESTING
    275     if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    276         fprintf(stderr, "Marking pixel image %d, pixel %d,%d for rejection\n", source, x, y);
    277     }
    278 #endif
     278    CHECKPIX(x, y, "Marking pixel image %d, pixel %d,%d for rejection\n", source, x, y);
    279279    pmStackData *data = inputs->data[source]; // Stack data of interest
    280280    if (!data) {
     
    290290static void combineExtract(int *num,                        // Number of good pixels
    291291                           bool *suspect,                   // Any suspect pixels?
     292                           psImageMaskType *badMask,        // OR of all bad (masked) pixels
     293                           psImageMaskType *goodMask,       // OR of all good (unmasked) pixels
    292294                           combineBuffer *buffer, // Buffer with vectors
    293295                           psImage *image, // Combined image, for output
     
    300302                           const psVector *reject, // Indices of pixels to reject, or NULL
    301303                           int x, int y, // Coordinates of interest; frame of output image
    302                            psImageMaskType maskVal, // Value to mask
    303                            psImageMaskType maskSuspect // Value to suspect
     304                           psImageMaskType badMaskBits, // Value to mask as 'bad'
     305                           psImageMaskType suspectMaskBits // Value to mask as 'suspect'
    304306                           )
    305307{
     
    322324    }
    323325
     326    // mask values to store possible mask combinations
     327    *badMask = 0;
     328    *goodMask = 0xffff;
     329
    324330    // Extract the pixel and mask data
    325331    int numGood = 0;                    // Number of good pixels
     
    328334        // should be because of how pixelMapGenerate works
    329335        if (reject && reject->data.U16[j] == i) {
     336            // pixels can be rejected because:
     337            // 1) only 1 input pixel (and 'safe' is true)
     338            // 2) only 2 input pixels were available and they were mutually inconsistent (or variance info was missing)
     339            // 3) NOTE : ifdef'ed out code for 3 inputs case
     340            // 4) outlier from sample of N pixels
     341            // 5) some of these may have been suspect.
     342            // XXX raise a specific mask bit for these (currently results in BLANK)
    330343            j++;
     344
     345            pmStackData *data = inputs->data[i]; // Stack data of interest
     346            if (data) {
     347                int xIn = x - data->readout->col0, yIn = y - data->readout->row0; // Coordinates on input readout
     348                psImage *mask = data->readout->mask; // Mask of interest
     349                *badMask |= mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn]; // save the bad bits used
     350                CHECKPIX(x, y, "reject: adding bit to mask: %d : %x (badMask = %x)\n", i, mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn], *badMask);
     351            } else {
     352                CHECKPIX(x, y, "reject: no item in data (badMask = %x)\n", *badMask);
     353            }
    331354            continue;
    332355        }
     
    334357        pmStackData *data = inputs->data[i]; // Stack data of interest
    335358        if (!data) {
     359            CHECKPIX(x, y, "skip: no item in data (badMask = %x)\n", *badMask);
    336360            continue;
    337361        }
     
    339363        int xIn = x - data->readout->col0, yIn = y - data->readout->row0; // Coordinates on input readout
    340364        psImage *mask = data->readout->mask; // Mask of interest
    341         if (mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & maskVal) {
     365        if (mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & badMaskBits) {
     366            *badMask |= (mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & badMaskBits); // save the bad bits used
     367            CHECKPIX(x, y, "skip: adding bit to mask: %d : %x (badMask = %x)\n", i, mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn], *badMask);
    342368            continue;
    343369        }
    344370
    345         pixelSuspects->data.U8[numGood] = mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & maskSuspect ?
    346             true : false;
     371        pixelSuspects->data.U8[numGood] = mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & suspectMaskBits ? true : false;
     372
     373        *goodMask &= mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn]; // save the mask bits still used
    347374
    348375        psImage *image = data->readout->image; // Image of interest
     
    356383        pixelSources->data.U16[numGood] = i;
    357384        numGood++;
     385
     386        CHECKPIX(x, y, "keep: %d : %x (badMask = %x)\n", i, mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn], *badMask);
    358387    }
    359388    pixelData->n = numGood;
     
    370399    if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    371400        for (int i = 0; i < numGood; i++) {
    372             fprintf(stderr, "Input %d, pixel %d,%d (%" PRIu16 "): %f %f (%f) %f %f %d\n",
     401            fprintf(stderr, "Input %d, pixel %d,%d (%" PRIu16 "): %f %f (%f) %f %f %d %x %x -> %x %x\n",
    373402                    i, x, y, pixelSources->data.U16[i], pixelData->data.F32[i], pixelVariances->data.F32[i],
    374403                    addVariance->data.F32[i], pixelWeights->data.F32[i], pixelExps->data.F32[i],
    375                     pixelSuspects->data.U8[i]);
     404                    pixelSuspects->data.U8[i], badMaskBits, suspectMaskBits, *badMask, *goodMask);
    376405        }
    377406    }
     
    380409    return;
    381410}
    382 
    383411
    384412// Combine pixels
     
    392420                          combineBuffer *buffer, // Buffer with vectors
    393421                          int x, int y, // Coordinates of interest; frame of output image
    394                           psImageMaskType bad, // Value for bad pixels
     422                          psImageMaskType blankMask, // Value for empty pixels
     423                          psImageMaskType badMask, // Value for bad pixels
     424                          psImageMaskType goodMask, // Value for good pixels
    395425                          bool safe,           // Safe combination?
    396426                          float invTotalWeight    // Inverse of total weight for all inputs
     
    404434    // Default option is that the pixel is bad
    405435    float imageValue = NAN, varianceValue = NAN; // Value for combined image and variance map
    406     psImageMaskType maskValue = bad;    // Value for combined mask
    407436    float expValue = 0.0, expWeightValue = NAN; // Exposure value (straight, and weighted)
     437
     438    // default output mask value.  badMask is the OR of all unused input pixels. 
     439    // if there are no input pixels, this value will be 0, in which case we want to set the output pixel to BLANK.
     440    // if there are only good input pixels, and they do not result in a valid pixel, we still want to set this to BLANK.
     441    psImageMaskType maskValue = badMask ? badMask : blankMask;    // Value for combined mask
     442
     443    CHECKPIX(x, y, "bad vs good : %x %x %x\n", maskValue, badMask, blankMask);
    408444
    409445    switch (num) {
    410446      case 0: {
    411447          // Nothing to combine: it's bad
    412 #ifdef TESTING
    413           if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    414               fprintf(stderr, "No inputs to combine, pixel %d,%d is bad.\n", x, y);
    415           }
    416 #endif
     448          CHECKPIX(x, y, "No inputs to combine, pixel %d,%d is bad.\n", x, y);
    417449          break;
    418450      }
     
    428460                  expWeightValue = pixelExps->data.F32[0];
    429461              }
    430               maskValue = 0;
    431 #ifdef TESTING
    432               if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    433                   fprintf(stderr, "Single input to combine, safety off, pixel %d,%d --> %f\n",
    434                           x, y, imageValue);
    435               }
    436 #endif
    437           }
    438 #ifdef TESTING
    439           else if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    440               fprintf(stderr, "Single input to combine, safety on, pixel %d,%d is bad.\n", x, y);
    441           }
    442 #endif
     462              maskValue = goodMask;
     463              CHECKPIX(x, y, "Single input to combine, safety off, pixel %d,%d --> %f\n", x, y, imageValue);
     464          } else {
     465              CHECKPIX(x, y, "Single input to combine, safety on, pixel %d,%d is bad.\n", x, y);
     466          }
    443467          break;
    444468      }
     
    446470          // Automatically accept the mean of the pixels only if we're not playing safe
    447471          if (!safe) {
    448               if (combinationMeanVariance(&imageValue, &varianceValue, &expValue, &expWeightValue,
    449                                           pixelData, pixelVariances, pixelExps, pixelWeights)) {
    450                   maskValue = 0;
    451 #ifdef TESTING
    452                   if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    453                       fprintf(stderr, "Two inputs to combine using unsafe, pixel %d,%d --> %f %f\n",
    454                               x, y, imageValue, varianceValue);
    455                   }
    456 #endif
     472              if (combinationMeanVariance(&imageValue, &varianceValue, &expValue, &expWeightValue, pixelData, pixelVariances, pixelExps, pixelWeights)) {
     473                  maskValue = goodMask;
     474                  CHECKPIX(x, y, "Two inputs to combine using unsafe, pixel %d,%d --> %f %f\n", x, y, imageValue, varianceValue);
    457475              }
     476          } else {
     477              CHECKPIX(x, y, "Two inputs to combine, safety on, pixel %d,%d is bad\n", x, y);
    458478          }
    459 #ifdef TESTING
    460           else {
    461               if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    462                   fprintf(stderr, "Two inputs to combine, safety on, pixel %d,%d is bad\n", x, y);
    463               }
    464           }
    465 #endif
    466479          break;
    467480      }
     
    472485              break;
    473486          }
    474           maskValue = 0;
    475 #ifdef TESTING
    476           if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    477               fprintf(stderr, "Combined inputs, pixel %d,%d --> %f %f\n", x, y, imageValue, varianceValue);
    478           }
    479 #endif
     487          maskValue = goodMask;
     488          CHECKPIX(x, y, "Combined inputs, pixel %d,%d --> %f %f\n", x, y, imageValue, varianceValue);
    480489          break;
    481490      }
     
    522531    int numIter = PS_MAX(iter * num, 1); // Number of iterations
    523532
    524 #ifdef TESTING
    525     if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    526         fprintf(stderr, "Testing pixel %d,%d: %d %f %f %f %d %d\n",
    527                 x, y, numIter, rej, sys, olympic, useVariance, safe);
    528     }
    529 #endif
     533    CHECKPIX(x, y, "Testing pixel %d,%d: %d %f %f %f %d %d\n", x, y, numIter, rej, sys, olympic, useVariance, safe);
    530534
    531535    psVector *pixelData = buffer->pixels; // Values for the pixel of interest
     
    548552            // Systematic error contributes to the rejection level
    549553            float sysVar = PS_SQR(sys * pixelData->data.F32[i]);
    550 #ifdef TESTING
    551             // Correct variance for comparison against weighted mean including itself
    552             float compare = 1.0 - pixelWeights->data.F32[i] / sumWeights;
    553             if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    554                 fprintf(stderr, "Variance %d (%d), pixel %d,%d: %f %f %f\n", i, pixelSources->data.U16[i],
    555                         x, y, pixelVariances->data.F32[i], sysVar, compare);
    556             }
    557 #endif
     554            CHECKPIX(x, y, "Variance %d (%d), pixel %d,%d: %f %f %f\n",
     555                     i, pixelSources->data.U16[i], x, y,
     556                     pixelVariances->data.F32[i], sysVar, 1.0 - pixelWeights->data.F32[i] / sumWeights);
    558557            pixelLimits->data.F32[i] = rej2 * (pixelVariances->data.F32[i] + sysVar);
    559558        }
     
    593592                          combineMarkInspect(inputs, x, y, pixelSources->data.U16[1]);
    594593                      }
    595 #ifdef TESTING
    596                       if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    597                           fprintf(stderr, "Flagged both inputs for pixel %d,%d (%f > %f x %f\n)",
    598                                   x, y, diff, rej, sqrtf(sigma2));
    599                       }
    600 #endif
     594                      CHECKPIX(x, y, "Flagged both inputs for pixel %d,%d (%f > %f x %f\n)", x, y, diff, rej, sqrtf(sigma2));
    601595                  }
    602596              } else if (i == 0 && safe) {
     
    708702                  float median = combinationWeightedOlympic(pixelData, pixelWeights,
    709703                                                            olympic, buffer->sort); // Median for stack
    710 #ifdef TESTING
    711                   if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    712                       fprintf(stderr, "Flag with variance pixel %d,%d: median = %f\n", x, y, median);
    713                   }
    714 #endif
     704                  CHECKPIX(x, y, "Flag with variance pixel %d,%d: median = %f\n", x, y, median);
    715705                  float worst = -INFINITY; // Largest deviation
    716706                  for (int j = 0; j < num; j++) {
    717707                      float diff = pixelData->data.F32[j] - median; // Difference from expected
    718 #ifdef TESTING
    719                       if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    720                           fprintf(stderr, "Testing input %d for pixel %d,%d: %f\n", j, x, y, diff);
    721                       }
    722 #endif
     708                      CHECKPIX(x, y, "Testing input %d for pixel %d,%d: %f\n", j, x, y, diff);
    723709
    724710                      // Comparing squares --- cheaper than lots of sqrts
     
    737723                  combinationMedianStdev(&median, &stdev, pixelData, buffer->sort);
    738724                  float limit = rej * stdev; // Rejection limit
    739 #ifdef TESTING
    740                   if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    741                       fprintf(stderr,
    742                               "Flag without variance pixel %d,%d; median = %f, stdev = %f, limit = %f\n",
    743                               x, y, median, stdev, limit);
    744                   }
    745 #endif
     725                  CHECKPIX(x, y, "Flag without variance pixel %d,%d; median = %f, stdev = %f, limit = %f\n", x, y, median, stdev, limit);
    746726                  float worst = -INFINITY; // Largest deviation
    747727                  for (int j = 0; j < num; j++) {
     
    763743        if (maskIndex >= 0) {
    764744            if (suspect) {
    765 #ifdef TESTING
    766                 if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    767                     fprintf(stderr, "Throwing out all suspect pixels for %d,%d\n", x, y);
    768                 }
    769 #endif
     745                CHECKPIX(x, y, "Throwing out all suspect pixels for %d,%d\n", x, y);
    770746                // Throw out all suspect pixels
    771747                int numGood = 0;        // Number of good pixels
     
    796772            } else {
    797773                // Throw out masked pixel
    798 #ifdef TESTING
    799                 if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    800                     fprintf(stderr, "Throwing out input %d for pixel %d,%d\n", maskIndex, x, y);
    801                 }
    802 #endif
     774                CHECKPIX(x, y, "Throwing out input %d for pixel %d,%d\n", maskIndex, x, y);
    803775                combineMarkInspect(inputs, x, y, pixelSources->data.U16[maskIndex]);
    804776                int numGood = 0;        // Number of good pixels
     
    999971
    1000972/// Stack input images
    1001 bool pmStackCombine(pmReadout *combined, pmReadout *expmaps, psArray *input,
    1002                     psImageMaskType maskVal, psImageMaskType maskSuspect,
    1003                     psImageMaskType bad, int kernelSize,
    1004                     float iter, float rej, float sys, float olympic,
    1005                     bool useVariance, bool safe, bool rejection)
     973bool pmStackCombine(
     974    pmReadout *combined,                // output stacked readout
     975    pmReadout *expmaps,                 // output exposure map information
     976    psArray *input,                     // input exposures
     977    psImageMaskType badMaskBits,        // treat these bits as 'bad'
     978    psImageMaskType suspectMaskBits,    // treat these bits as 'suspect'
     979    psImageMaskType blankMaskBits,      // use this mask value for pixels missing input data (distinguish between Ninput = 0 and Ngood = 0?)
     980    int kernelSize,
     981    float iter,
     982    float rej,
     983    float sys,
     984    float olympic,
     985    bool useVariance,
     986    bool safe,
     987    bool rejection)
    1006988{
    1007989    bool haveVariances;                 // Do we have the variance maps?
     
    1012994    }
    1013995    PS_ASSERT_INT_NONNEGATIVE(kernelSize, false);
    1014     PS_ASSERT_INT_POSITIVE(bad, false);
     996    PS_ASSERT_INT_POSITIVE(blankMaskBits, false);
    1015997    if (isnan(rej)) {
    1016998        PS_ASSERT_FLOAT_EQUAL(iter, 0, false);
     
    10701052    }
    10711053    psFree(stack);
    1072     pmReadoutUpdateSize(combined, minInputCols, minInputRows, xSize, ySize, true, true, bad);
     1054    pmReadoutUpdateSize(combined, minInputCols, minInputRows, xSize, ySize, true, true, blankMaskBits);
    10731055    psTrace("psModules.imcombine", 1, "Have for combination [%d:%d,%d:%d] (%dx%d)\n",
    10741056            minInputCols, maxInputCols, minInputRows, maxInputRows, xSize, ySize);
     
    11311113    for (int y = minInputRows; y < maxInputRows; y++) {
    11321114        for (int x = minInputCols; x < maxInputCols; x++) {
     1115            CHECKPIX(x, y, "Combining pixel %d,%d: %x %x %f %f %f %f %d %d %d\n", x, y, badMaskBits, blankMaskBits, iter, rej, sys, olympic, useVariance, safe, rejection);
     1116
    11331117#ifdef TESTING
    1134             if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    1135                 fprintf(stderr, "Combining pixel %d,%d: %x %x %f %f %f %f %d %d %d\n",
    1136                         x, y, maskVal, bad, iter, rej, sys, olympic, useVariance, safe, rejection);
    1137             }
    1138 #endif
     1118            if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
     1119                fprintf (stderr, "combine\n");
     1120            }
     1121# endif
     1122
    11391123            psVector *reject = NULL; // Images to reject for this pixel
    11401124            if (rejection) {
     
    11541138#endif
    11551139            }
    1156 
    1157             int num;                    // Number of good pixels
    1158             bool suspect;               // Suspect pixels in stack?
    1159             combineExtract(&num, &suspect, buffer, combinedImage, combinedMask, combinedVariance,
    1160                            input, weights, exps, addVariance, reject, x, y, maskVal, maskSuspect);
    1161             combinePixels(combinedImage, combinedMask, combinedVariance, exp, expnum, expweight,
    1162                           num, buffer, x, y, bad, safe, totalExpWeight);
     1140 
     1141            int num;                      // Number of good pixels
     1142            bool suspect;                 // Suspect pixels in stack?
     1143            psImageMaskType badMask = 0;  // OR of mask bits in all bad input pixels
     1144            psImageMaskType goodMask = 0; // OR of mask bits in all good input pixels
     1145            combineExtract(&num, &suspect, &badMask, &goodMask, buffer, combinedImage, combinedMask, combinedVariance, input, weights, exps, addVariance, reject, x, y, badMaskBits, suspectMaskBits);
     1146            combinePixels(combinedImage, combinedMask, combinedVariance, exp, expnum, expweight, num, buffer, x, y, blankMaskBits, badMask, goodMask, safe, totalExpWeight);
    11631147
    11641148            if (iter > 0) {
  • branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionMatch.c

    r29170 r29215  
    275275}
    276276
    277 bool pmSubtractionMaskInvalid (const pmReadout *readout, psImageMaskType maskVal) {
    278 
    279     if (!readout) return true;
    280 
    281     psImage *image = readout->image;
    282     psImage *mask  = readout->mask;
    283     psImage *variance = readout->variance;
    284     for (int y = 0; y < image->numRows; y++) {
    285         for (int x = 0; x < image->numCols; x++) {
    286             if (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskVal) continue;
    287             bool valid = false;
    288             valid = isfinite(image->data.F32[y][x]);
    289             if (variance) {
    290                 valid &= isfinite(variance->data.F32[y][x]);
    291             }
    292             if (valid) continue;
    293             mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] = maskVal;
    294         }
    295     }
    296 
    297     return true;
    298 }
     277// bool pmSubtractionMaskInvalid (const pmReadout *readout, psImageMaskType maskVal) {
     278//
     279//     if (!readout) return true;
     280//
     281//     psImage *image = readout->image;
     282//     psImage *mask  = readout->mask;
     283//     psImage *variance = readout->variance;
     284//     for (int y = 0; y < image->numRows; y++) {
     285//         for (int x = 0; x < image->numCols; x++) {
     286//             if (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskVal) continue;
     287//             bool valid = false;
     288//             valid = isfinite(image->data.F32[y][x]);
     289//             if (variance) {
     290//                 valid &= isfinite(variance->data.F32[y][x]);
     291//             }
     292//             if (valid) continue;
     293//             mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] = maskVal;
     294//         }
     295//     }
     296//
     297//     return true;
     298// }
    299299
    300300bool pmSubtractionMatchPrecalc(pmReadout *conv1, pmReadout *conv2, const pmReadout *ro1, const pmReadout *ro2,
     
    387387    }
    388388
    389     pmSubtractionMaskInvalid(ro1, maskVal);
    390     pmSubtractionMaskInvalid(ro2, maskVal);
     389    // XXX this is done before calling this function
     390    // pmSubtractionMaskInvalid(ro1, maskVal);
     391    // pmSubtractionMaskInvalid(ro2, maskVal);
    391392
    392393    // General background subtraction, since this is done in pmSubtractionMatch
     
    560561    memCheck("start");
    561562
    562     pmSubtractionMaskInvalid(ro1, maskVal);
    563     pmSubtractionMaskInvalid(ro2, maskVal);
     563    // pmSubtractionMaskInvalid(ro1, maskVal);
     564    // pmSubtractionMaskInvalid(ro2, maskVal);
    564565
    565566    psRegion bounds = psRegionSet(NAN, NAN, NAN, NAN); // Bounds of valid pixels
  • branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionStamps.c

    r29170 r29215  
    106106
    107107    if (xMax < xMin) {
    108       fprintf (stderr, "%f,%f : x-border\n", xRaw, yRaw);
    109       return false;
     108        // fprintf (stderr, "%f,%f : x-border\n", xRaw, yRaw);
     109        return false;
    110110    }
    111111    if (yMax < yMin) {
    112       fprintf (stderr, "%f,%f : y-border\n", xRaw, yRaw);
    113       return false;
     112        // fprintf (stderr, "%f,%f : y-border\n", xRaw, yRaw);
     113        return false;
    114114    }
    115115
     
    458458                // Take stamps off the top of the (sorted) list
    459459                for (int j = xList->n - 1; j >= 0 && !goodStamp; j--) {
    460                   fprintf (stderr, "%d : xList: %ld elements\n", i, xList->n);
     460                    // fprintf (stderr, "%d : xList: %ld elements\n", i, xList->n);
    461461                    // Chop off the top of the list
    462462                    xList->n = j;
     
    706706
    707707    // storage vector for flux data
    708     psStats *stats = psStatsAlloc (PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV);
     708    psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
    709709    psVector *flux1 = psVectorAllocEmpty(2*stamps->num, PS_TYPE_F32);
    710710    psVector *flux2 = psVectorAllocEmpty(2*stamps->num, PS_TYPE_F32);
     
    734734                psAbort ("failed to generate stats");
    735735            }
    736             float f1 = stats->robustMedian;
     736            float f1 = stats->sampleMedian;
    737737
    738738            psStatsInit (stats);
     
    740740                psAbort ("failed to generate stats");
    741741            }
    742             float f2 = stats->robustMedian;
     742            float f2 = stats->sampleMedian;
    743743
    744744            stamps->window1->kernel[y][x] = f1;
     
    755755        }
    756756    }
     757
     758#if 0
     759    {
     760        psFits *fits = NULL;
     761        fits = psFitsOpen ("window1.fits", "w");
     762        psFitsWriteImage (fits, NULL, stamps->window1->image, 0, NULL);
     763        psFitsClose (fits);
     764        fits = psFitsOpen ("window2.fits", "w");
     765        psFitsWriteImage (fits, NULL, stamps->window2->image, 0, NULL);
     766        psFitsClose (fits);
     767    }
     768#endif
    757769
    758770    psTrace("psModules.imcombine", 3, "Window total (1): %f, threshold: %f\n", sum1, (1.0 - stamps->normFrac) * sum1);
Note: See TracChangeset for help on using the changeset viewer.