IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 13390


Ignore:
Timestamp:
May 15, 2007, 3:04:43 PM (19 years ago)
Author:
Paul Price
Message:

Spatial variation working for both POIS and ISIS kernels.

Location:
trunk/psModules/src/imcombine
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/psModules/src/imcombine/pmSubtraction.c

    r13389 r13390  
    44 *  @author GLG, MHPCC
    55 *
    6  *  @version $Revision: 1.7 $ $Name: not supported by cvs2svn $
    7  *  @date $Date: 2007-05-16 00:01:03 $
     6 *  @version $Revision: 1.8 $ $Name: not supported by cvs2svn $
     7 *  @date $Date: 2007-05-16 01:04:43 $
    88 *
    99 *  Copyright 2004-2007 Institute for Astronomy, University of Hawaii
     
    5858                              const pmSubtractionKernels *kernels, // Kernel basis functions
    5959                              const psImage *polyValues, // Spatial polynomial values
    60                               float (*weightFunc)(float value) // Function for weighting
     60                              double (*weightFunc)(double value) // Function for weighting
    6161    )
    6262{
     
    8888              int v = kernels->v->data.S32[i]; // Offset in y
    8989              kernel->kernel[v][u] += weightFunc(polyValue * solution->data.F64[i]);
    90               if (i != kernels->subIndex) {
     90              if (kernels->spatialOrder > 0 && i != kernels->subIndex) {
    9191                  // The (0,0) element is subtracted from most kernels to preserve photometric scaling
    92                   kernel->kernel[0][0] -= weightFunc(solution->data.F64[i]);
     92                  kernel->kernel[0][0] += weightFunc(-solution->data.F64[i]);
    9393              }
    9494              break;
     
    9696          case PM_SUBTRACTION_KERNEL_ISIS: {
    9797              psKernel *preCalc = kernels->preCalc->data[i];// Precalculated values
     98              psKernel *subKernel = kernels->preCalc->data[kernels->subIndex]; // Kernel to subtract
    9899              for (int v = -size; v <= size; v++) {
    99100                  for (int u = -size; u <= size; u++) {
    100                       kernel->kernel[v][u] += weightFunc(polyValue * preCalc->kernel[v][u]);
     101                      kernel->kernel[v][u] += weightFunc(solution->data.F64[i] * polyValue *
     102                                                         preCalc->kernel[v][u]);
     103                      // The (0,0) kernel is subtracted from other kernels to preserve photometric scaling
     104                      if (kernels->spatialOrder > 0 && i != kernels->subIndex) {
     105                          kernel->kernel[v][u] += weightFunc(-solution->data.F64[i] *
     106                                                             subKernel->kernel[v][u]);
     107                      }
    101108                  }
    102109              }
     
    117124                                   const psImage *image, // Image to convolve
    118125                                   const psImage *polyValues, // Spatial polynomial values
    119                                    float (*weightFunc)(float value) // Function for weighting
     126                                   double (*weightFunc)(double value) // Function for weighting
    120127    )
    121128{
     
    130137          int v = kernels->v->data.S32[index]; // Offset in y
    131138          double value = weightFunc(polyValue) * image->data.F32[y + v][x + u]; // Value of convolution
    132           if (index != kernels->subIndex) {
     139          if (kernels->spatialOrder > 0 && index != kernels->subIndex) {
    133140              // The (0,0) element is subtracted from most kernels to preserve photometric scaling
    134               value -= image->data.F32[y][x];
     141              value += weightFunc(-1.0) * image->data.F32[y][x];
    135142          }
    136143          return value;
     
    138145      case PM_SUBTRACTION_KERNEL_ISIS: {
    139146          psKernel *kernel = kernels->preCalc->data[index]; // The convolution kernel
     147          psKernel *subKernel = kernels->preCalc->data[kernels->subIndex]; // Kernel to subtract
    140148          int size = kernels->size;     // Kernel half-size
    141149          double sum = 0.0;             // Accumulated sum from convolution
     150          double sub = 0.0;             // Accumulated sum to subtract
    142151          for (int v = -size; v <= size; v++) {
    143152              for (int u = -size; u <= size; u++) {
    144153                  sum += weightFunc(kernel->kernel[v][u]) * image->data.F32[y + v][x + u];
     154                  // The (0,0) kernel is subtracted from other kernels to preserve photometric scaling
     155                  if (kernels->spatialOrder > 0 && index != kernels->subIndex) {
     156                      sub += weightFunc(subKernel->kernel[v][u]) * image->data.F32[y + v][x + u];
     157                  }
    145158              }
    146159          }
    147           return weightFunc(polyValue) * sum;
     160          return weightFunc(polyValue) * sum + weightFunc(-1.0) * sub;
    148161      }
    149162      default:
     
    154167
    155168// Weighting function for use with convolvePixel: no weighting applied, suitable for combining image pixels
    156 static inline float imageWeighting(float value)
     169static inline double imageWeighting(double value)
    157170{
    158171    return value;
     
    160173
    161174// Weighting function for use with convolvePixel: weighting suitable for combining variances
    162 static inline float varianceWeighting(float value)
     175static inline double varianceWeighting(double value)
    163176{
    164177    return PS_SQR(value);
  • trunk/psModules/src/imcombine/pmSubtractionKernels.c

    r13389 r13390  
    137137                }
    138138
     139#if 0
     140                // Subtract a particular kernel in order to preserve photometric calibration across image
     141                if (spatialOrder > 0 && index != kernels->subIndex) {
     142                    psKernel *subKernel = kernels->preCalc->data[kernels->subIndex]; // Kernel to subtract
     143                    (void)psBinaryOp(preCalc->image, preCalc->image, "-", subKernel->image);
     144                }
     145#endif
     146
    139147                // Iterate over spatial order.  This loop creates the terms for
    140148                // x^xOrder * y^yOrder  such that (xOrder+yOrder) <= spatialOrder.
     
    153161                }
    154162
    155 
    156                 if (psTraceGetLevel("psModules.imcombine.kernel") >= 10) {
    157                     psString kernelName = NULL;
    158                     psStringAppend(&kernelName, "kernel%d.fits", index);
    159                     psFits *kernelFile = psFitsOpen(kernelName, "w");
    160                     psFree(kernelName);
    161                     psFitsWriteImage(kernelFile, NULL, preCalc->image, 0, NULL);
    162                     psFitsClose(kernelFile);
    163                 }
    164 
    165163                psFree(preCalc);        // Drop reference
    166164            }
     
    168166    }
    169167
    170     // Subtract a particular kernel in order to preserve photometric calibration across image
    171168    kernels->subIndex = 0;
    172169    assert(kernels->u->data.S32[kernels->subIndex] == 0 &&
     
    175172           kernels->yOrder->data.S32[kernels->subIndex] == 0);
    176173
    177     psKernel *subKernel = kernels->preCalc->data[kernels->subIndex]; // Kernel to subtract
    178     for (int i = 0; i < kernels->subIndex; i++) {
    179         psKernel *kernel = kernels->preCalc->data[i]; // Kernel of interest
    180         (void)psBinaryOp(kernel->image, kernel->image, "-", subKernel->image);
    181     }
    182     for (int i = kernels->subIndex + 1; i < num; i++) {
    183         psKernel *kernel = kernels->preCalc->data[i]; // Kernel of interest
    184         (void)psBinaryOp(kernel->image, kernel->image, "-", subKernel->image);
    185     }
    186 
     174    if (psTraceGetLevel("psModules.imcombine.kernel") >= 10) {
     175        for (int i = 0; i < num; i++) {
     176            psKernel *kernel = kernels->preCalc->data[i]; // Kernel of interest
     177            psString kernelName = NULL;
     178            psStringAppend(&kernelName, "kernel%03d.fits", i);
     179            psFits *kernelFile = psFitsOpen(kernelName, "w");
     180            psFree(kernelName);
     181            psFitsWriteImage(kernelFile, NULL, kernel->image, 0, NULL);
     182            psFitsClose(kernelFile);
     183        }
     184    }
    187185
    188186    return kernels;
Note: See TracChangeset for help on using the changeset viewer.