Changeset 6761

Timestamp:

Apr 4, 2006, 8:41:38 AM (20 years ago)

Author:

magnier

Message:

finished first cut at fringe stats code

Location:

branches/rel10_ifa/psModules/src/detrend

Files:

• pmFringeStats.c (modified) (3 diffs)
• pmFringeStats.h (modified) (2 diffs)

branches/rel10_ifa/psModules/src/detrend/pmFringeStats.c

@@ -3 +3 @@
  *  @author Eugene Magnier, IfA
  *
- *  @version $Revision: 1.1.2.1 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2006-04-03 06:12:10 $
+ *  @version $Revision: 1.1.2.2 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-04 18:41:38 $
  *
  *  Copyright 2004 IfA
@@ -44 +44 @@
 }
 
+static pmFringeScaleFree (pmFringeScale *scale)
+{
+    psFree (scale->coeff);
+    psFree (scale->coeffErr);
+    return;
+}
+
+pmFringeScale *pmFringeScaleAlloc (int nFringeFrames)
+{
+    pmFringeScale *scale = psAlloc(sizeof(pmFringeScale));
+    (void)psMemSetDeallocator(scale, (psFreeFunc)pmFringeScaleFree);
+
+    scale->nFringeFrames = nFringeFrames;
+    scale->coeff = psVectorAlloc (nFringeFrames + 1, PS_TYPE_F32);
+    scale->coeffErr = psVectorAlloc (nFringeFrames + 1, PS_TYPE_F32);
+
+    return scale;
+}
+
 bool pmFringeStatsCreatePoints (pmFringeStats *fringe, psImage *image)
 {
@@ -117 +136 @@
     return true;
 }
-psImage *pmFringeCorrect(
-    psImage *out,
-    psMetadata *info,
-    psImage *science,
-    psArray *fringeImage,
-    psArray *fringeStats
-);
-
+
+// XXX include the fringe error (fringe->df) in the fit?
+pmFringeScale *pmFringeScaleMeasure (pmFringeStats *science, psArray *fringes)
+{
+
+    double sum;
+    psVector *v1;
+    psVector *v2;
+
+    pmFringeScale *scale = pmFringeScaleAlloc (fringes->n);
+
+    int nCof = fringes->n + 1;
+    int nPts = science->f->n;
+
+    psImage *A = psImageAlloc (nCof, nCof, PS_TYPE_F64);
+    psVector *B = psVectorAlloc (nCod, PS_TYPE_F64);
+
+    for (int i = 0; i < nCof; i++) {
+        if (i == 0) {
+            p1 = NULL;
+        } else {
+            v1 = fringes->data[i - 1];
+            p1 = v1->data.F32;
+        }
+        for (int j = 0; j < nCof; j++) {
+            if (j == 0) {
+                p2 = NULL;
+            } else {
+                v2 = fringes->data[j - 1];
+                p2 = v2->data.F32;
+            }
+
+            sum = 0;
+            for (int I = 0; I < nPts; I++) {
+                f1 = (p1 == NULL) ? 1.0 : p1[I];
+                f2 = (p2 == NULL) ? 1.0 : p2[I];
+                sum += f1*f2;
+            }
+            A->data.F32[i][j] = sum;
+        }
+
+        sum = 0;
+        p2 = science->f->data.F32;
+        for (int I = 0; I < nPts; I++) {
+            f1 = (p1 == NULL) ? 1.0 : p1[I];
+            f2 = p2[I];
+            sum += f1*f2;
+        }
+        B->data.F32[i] = sum;
+    }
+
+    if (! psGaussJordan(A, B)) {
+        psError(PS_ERR_UNKNOWN, false, "Could not solve linear equations.  Returning NULL.\n");
+    }
+
+
+    for (int i = 0; i < nCof; i++) {
+        scale->coeff->data.F32[i] = B->data.F32[i];
+        scale->coeffErr->data.F32[i] = sqrt(A->data.F32[i][i]);
+    }
+
+    return scale;
+}
+
+// XXX note that this modifies the input fringe images
+psImage *pmFringeCorrect(psImage *out, psMetadata *info, psImage *science, psArray *fringeImage, psArray *fringeStats)
+{
+
+    // measure the fringe stats
+    pmFringeScale *scale = pmFringeScaleMeasure (science, fringeStats);
+
+    if (out != NULL) {
+        psImageInit (out, 0.0);
+    }
+
+    // build the fringe correction image
+    // XXX we could save data space by making the first image the output image
+    for (int i = 0; i < fringeImage->n; i++) {
+
+        // rescale the fringe image
+        psBinaryOp (fringeImage->data[i], fringeImage->data[i], "*", psScalarAlloc (scale->coeff->data.F32[i+1]));
+
+        // sum together
+        out = psBinaryOp (out, out, "+", fringeImage->data[i]);
+    }
+
+    // subtract the resulting fringe frame
+    out = psBinaryOp (out, science, "-", out);
+
+    return out;
+}
+

branches/rel10_ifa/psModules/src/detrend/pmFringeStats.h

@@ -5 +5 @@
  *  @author Eugene Magnier, IfA
  *
- *  @version $Revision: 1.1.2.1 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2006-04-03 06:12:10 $
+ *  @version $Revision: 1.1.2.2 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-04 18:41:38 $
  *
  *  Copyright 2004 IfA, University of Hawaii
@@ -41 +41 @@
 );
 
+/** the pmFringeScale structure defines the relationship between two fringe measurements
+ */
+typedef struct
+{
+    int nFringeFrames;
+    psVector *coeff;
+    psVector *coeffErr;
+}
+pmFringeScale;
+
 /** Measure the fringe stats for an image
  *