Changeset 17702 for branches/eam_branch_20080511/ppSim/src/ppSimMakeSky.c

Timestamp:

May 15, 2008, 12:01:36 PM (18 years ago)

Author:

eugene

Message:

working on getting the input image options right

File:

• branches/eam_branch_20080511/ppSim/src/ppSimMakeSky.c (modified) (4 diffs)

branches/eam_branch_20080511/ppSim/src/ppSimMakeSky.c

@@ -1 +1 @@
 # include "ppSim.h"
 
-// XXX add bounds to the inputs?
+// this function sets the skyRate to a value for the night sky (SKY.RATE or SKY.MAGS) or for a
+// flat-field image (FLAT.RATE).  Include a shutter correction and a scattered light source
+
 bool ppSimMakeSky (pmReadout *readout, psImage *expCorr, ppSimType type, pmConfig *config) {
 
@@ -16 +18 @@
 
     bool sky = psMetadataLookupBool(&status, recipe, "SKY"); // Generate a SKY flux?
-    if (!sky) return true;
  
-    float expTime     = psMetadataLookupF32(&status, recipe, "EXPTIME"); // Exposure time
-    float flatSigma   = psMetadataLookupF32(&status, recipe, "FLAT.SIGMA"); // Flat-field coefficient
-    float flatRate    = psMetadataLookupF32(&status, recipe, "FLAT.RATE"); // Flat-field rate
-    float shutterTime = psMetadataLookupF32(&status, recipe, "SHUTTER.TIME"); // Shutter time
-    float skyRate     = psMetadataLookupF32(&status, recipe, "SKY.RATE"); // Sky rate
+    float expTime      = psMetadataLookupF32(&status, recipe, "EXPTIME"); // Exposure time
+
+    float flatSigma    = psMetadataLookupF32(&status, recipe, "FLAT.SIGMA"); // Flat-field coefficient
+    float flatRate     = psMetadataLookupF32(&status, recipe, "FLAT.RATE"); // Flat-field rate
+    float shutterTime  = psMetadataLookupF32(&status, recipe, "SHUTTER.TIME"); // Shutter time
+    float scatterScale = psMetadataLookupF32(&status, recipe, "SCATTER.SCALE"); // scattered light term
+    float skyRate      = psMetadataLookupF32(&status, recipe, "SKY.RATE"); // Sky rate
     if (isnan(skyRate)) {
-        float zp      = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); assert (status);
-        float scale   = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE"); assert (status);
-        float skyMags = psMetadataLookupF32(&status, recipe, "SKY.MAGS");  assert (status);
+        float zp       = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); assert (status);
+        float scale    = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE"); assert (status);
+        float skyMags  = psMetadataLookupF32(&status, recipe, "SKY.MAGS");  assert (status);
         skyRate = scale * scale * ppSimMagToFlux (skyMags, zp);
+    }
+    if (type == PPSIM_TYPE_FLAT) {
+      skyRate = flatRate;
     }
 
@@ -44 +50 @@
     // Size of FPA
     psRegion *bounds = ppSimFPABounds (fpa);
+    int dXfpa = bounds->x1 - bounds->x0;
+    int dYfpa = bounds->y1 - bounds->y0;
 
     // Correct chip offsets so that boresight is in the middle of the FPA
-    x0Chip -= 0.5 * (bounds->x1 - bounds->x0);
-    y0Chip -= 0.5 * (bounds->y1 - bounds->y0);
+    x0Chip -= 0.5 * dXfpa;
+    y0Chip -= 0.5 * dYfpa;
 
     for (int y = 0; y < signal->numRows; y++) {
@@ -62 +70 @@
 
             // Gaussian flat-field over the FPA
-            float flatValue = expf(-0.5 / PS_SQR(flatSigma) * (PS_SQR(yFPA) + PS_SQR(xFPA))) /
-                flatSigma / sqrtf(2.0 * M_PI);
+            float flatValue = expf(-0.5 / PS_SQR(flatSigma) * (PS_SQR(yFPA) + PS_SQR(xFPA))) / flatSigma / sqrtf(2.0 * M_PI);
 
-            if (type == PPSIM_TYPE_FLAT) {
-                float flatFlux = flatRate * flatValue * realExpTime; // Flux from flat-field
-                signal->data.F32[y][x] += flatFlux;
-                variance->data.F32[y][x] += flatFlux;
-                continue;
-            }
+            float scatterRate = 0.0;
 
-            expCorr->data.F32[y][x] = realExpTime / expTime;
+            if (sky) {
+              // add a scattered light term to the flat-field images (no
+              if (type == PPSIM_TYPE_FLAT) {
+                float xF = 2.0*(xFPA / dXfpa) - 1.0;
+                scatterRate = scatterScale * PS_SQR(xF);
+              }
 
-            // Sky background
-            float skyFlux = skyRate * flatValue * realExpTime; // Flux from sky
-            signal->data.F32[y][x] += skyFlux;
-            variance->data.F32[y][x] += skyFlux;
+              // Sky background
+              float skyFlux = (skyRate * flatValue + scatterRate) * realExpTime; // Flux from sky
+              signal->data.F32[y][x] += skyFlux;
+              variance->data.F32[y][x] += skyFlux;
+            }
+
+            // used later to modify the star and galaxy photometry 
+            if (expCorr) {
+              expCorr->data.F32[y][x] = realExpTime / expTime;
+            }
 
             // TO DO: Add fringes