Changeset 18011 for trunk/ppSim/src/ppSimMakeSky.c

Timestamp:

Jun 8, 2008, 4:03:31 PM (18 years ago)

Author:

eugene

Message:

update from changes on eam_branch_20080511 : adding photometry of fake sources, force photometry; major cleanups

File:

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

trunk/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) {
 
@@ -15 +17 @@
     psMetadata *recipe = psMetadataLookupMetadata(&status, config->recipes, PPSIM_RECIPE); // Recipe
 
-    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
-    if (isnan(skyRate)) {
-        float zp      = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); assert (status);
-        float scale   = psMetadataLookupF32(&status, recipe, "SCALE");     assert (status);
-        float skyMags = psMetadataLookupF32(&status, recipe, "SKY.MAGS");  assert (status);
+    bool sky  = psMetadataLookupBool(&status, recipe, "SKY"); // Generate a SKY flux?
+    bool flat = psMetadataLookupBool(&status, recipe, "FLAT"); // Apply flat-field term?
+ 
+    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 scatterFrac  = psMetadataLookupF32(&status, recipe, "SCATTER.FRAC"); // scattered light fraction (max)
+    float skyRate      = psMetadataLookupF32(&status, recipe, "SKY.RATE"); // Sky rate
+    float skyMags      = psMetadataLookupF32(&status, recipe, "SKY.MAGS");  assert (status);
+    if (!isnan(skyMags)) {
+        float zp       = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); assert (status);
+        float scale    = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE"); assert (status);
         skyRate = scale * scale * ppSimMagToFlux (skyMags, zp);
+    }
+    if (type == PPSIM_TYPE_FLAT) {
+      skyRate = flatRate;
     }
 
@@ -37 +47 @@
     int yParityCell   = psMetadataLookupS32(&status, cell->concepts, "CELL.YPARITY");
 
-    int binning = psMetadataLookupS32(NULL, config->arguments, "BINNING"); // Binning in x and y
+    int binning = psMetadataLookupS32(NULL, recipe, "BINNING"); // Binning in x and y
 
     // 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++) {
@@ -58 +70 @@
             float realExpTime = expTime + shutterTime * (xFPA + yFPA + 2.0) / 4.0;
 
-            // 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);
+            // Gaussian flat-field over the FPA with flatValue = 1.0 at the field center
+            float flatValue = 1.0;
+            if (flat) {
+                // we make the flat-field have a response of 1.0 at the field center (like a vignetting)
+                flatValue = expf(-0.5 / PS_SQR(flatSigma) * (PS_SQR(yFPA) + PS_SQR(xFPA)));
+            }
 
-            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 = scatterFrac * 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