Changeset 14724

Timestamp:

Sep 1, 2007, 4:42:55 PM (19 years ago)

Author:

magnier

Message:

working on interpolation

File:

• branches/eam_branch_20070830/psLib/src/imageops/psImagePixelInterpolate.c (modified) (4 diffs)

branches/eam_branch_20070830/psLib/src/imageops/psImagePixelInterpolate.c

@@ -11 +11 @@
  *  @author Eugene Magnier, IfA
  *
- *  @version $Revision: 1.1.2.2 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2007-09-02 02:03:58 $
+ *  @version $Revision: 1.1.2.3 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2007-09-02 02:42:55 $
  *
  *  Copyright 2007 Institute for Astronomy, University of Hawaii
@@ -34 +34 @@
 #include "psPolynomial.h"
 #include "psPolynomialUtils.h"
+#include "psMinimizePolyFit.h"
 #include "psImage.h"
 #include "psImageInterpolate.h"
@@ -49 +50 @@
                     if (jx + ix >= mask->numCols) { continue; } \
                     /* do not test self */ \
-                    if (!jx && jy) { continue; } \
-                    if (mask->data.PS_TYPE_MASK_DATA[iy][ix] & maskVal) { continue; } \
+                    if (!jx && !jy) { continue; } \
+                    if (mask->data.PS_TYPE_MASK_DATA[iy+jy][ix+jx] & maskVal) { continue; } \
                     nGood ++; \
                 } \
@@ -115 +116 @@
                   // fit a quadratic to the valid neighbor pixels
                   // apply the quadratic to get the poor pixel value
-                  psPolynomial2D *poly = psImageBicubeFit (image, mask, maskVal, ix, iy);
+                  psPolynomial2D *poly = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, 2, 2);
+                  poly->mask[2][2] = 1;
+                  poly->mask[2][1] = 1;
+                  poly->mask[1][2] = 1;
+                  psVector *f = psVectorAlloc (9, PS_TYPE_F32);
+                  psVector *df = psVectorAlloc (9, PS_TYPE_F32);
+                  psVector *x = psVectorAlloc (9, PS_TYPE_F32);
+                  psVector *y = psVectorAlloc (9, PS_TYPE_F32);
+
+                  int n = 0;
+                  for (int jy = -1; jy <= +1; jy++) {
+                      // skip invalid pixels (jx < 0, >= Nx), etc
+                      if (jy + iy < 0) { continue; }
+                      if (jy + iy >= image->numRows) { continue; }
+                      for (int jx = -1; jx <= +1; jx++) {
+                          // skip invalid pixels (jx < 0, >= Nx), etc
+                          if (jx + ix < 0) { continue; } 
+                          if (jx + ix >= mask->numCols) { continue; } 
+                          /* skip self */
+                          if (!jx && !jy) { continue; } 
+                          // skip masked pixels
+                          if (mask->data.PS_TYPE_MASK_DATA[iy+jy][ix+jx] & maskVal) { continue; } \
+                          f->data.F32[n] = image->data.F32[iy+jy][ix+jx];
+                          // df->data.F32[n] = weight->data.F32[jy][jx];
+                          x->data.F32[n] = jx;
+                          y->data.F32[n] = jy;
+                          n++;
+                      }
+                  }
+                  // XXX set vector lengths here
+                  x->n = n;
+                  y->n = n;
+                  f->n = n;
+                  // df->n = n;
+                  psVectorFitPolynomial2D (poly, mask, maskValue, f, df, x, y);
                   image->data.F32[iy][ix] = poly->coeff[0][0];
                   break; }
 
               case PS_IMAGE_INTERPOLATE_LL: {
-                  // fit a plane to the valid neighbor pixels
-                  // apply the plane to get the poor pixel value
-                  psPolynomial2D *poly = psImageBilinearFit (image, mask, maskVal, ix - 1, iy - 1);
-                  image->data.F32[iy][ix] = poly->coeff[0][0] + poly->coeff[1][0] + poly->coeff[0][1];
+                  // fit a plane to the 3 pixels at (0,1),(1,0),(1,1), extend to pixel at (0,0)
+                  image->data.F32[iy][ix] = image->data.F32[iy+1][ix+0] - image->data.F32[iy+0][ix+1] - image->data.F32[iy+1][ix+1];
                   break; }
 
               case PS_IMAGE_INTERPOLATE_LR: {
-                  // fit a plane to the valid neighbor pixels
-                  // apply the plane to get the poor pixel value
-                  psPolynomial2D *poly = psImageBilinearFit (image, mask, maskVal, ix - 0, iy - 1);
-                  image->data.F32[iy][ix] = poly->coeff[0][0] + poly->coeff[0][1];
+                  // fit a plane to the 3 pixels at (0,1),(-1,0),(-1,1), extend to pixel at (0,0)
+                  image->data.F32[iy][ix] = image->data.F32[iy+1][ix+0] - image->data.F32[iy+0][ix-1] - image->data.F32[iy+1][ix-1];
                   break; }
 
               case PS_IMAGE_INTERPOLATE_UL: {
-                  // fit a plane to the valid neighbor pixels
-                  // apply the plane to get the poor pixel value
-                  psPolynomial2D *poly = psImageBilinearFit (image, mask, maskVal, ix - 1, iy - 0);
-                  image->data.F32[iy][ix] = poly->coeff[0][0] + poly->coeff[1][0];
+                  // fit a plane to the 3 pixels at (0,-1),(1,0),(1,-1), extend to pixel at (0,0)
+                  image->data.F32[iy][ix] = image->data.F32[iy-1][ix+0] - image->data.F32[iy+0][ix+1] - image->data.F32[iy-1][ix+1];
                   break; }
 
               case PS_IMAGE_INTERPOLATE_UR: {
-                  // fit a plane to the valid neighbor pixels
-                  // apply the plane to get the poor pixel value
-                  psPolynomial2D *poly = psImageBilinearFit (image, mask, maskVal, ix - 0, iy - 0);
-                  image->data.F32[iy][ix] = poly->coeff[0][0];
+                  // fit a plane to the 3 pixels at (0,-1),(-1,0),(-1,-1), extend to pixel at (0,0)
+                  image->data.F32[iy][ix] = image->data.F32[iy-1][ix+0] - image->data.F32[iy+0][ix-1] - image->data.F32[iy-1][ix-1];
                   break; }