Changeset 23335 for branches/eam_branches/eam_branch_20090312/psastro/src/psastroExtractAnalysis.c

Timestamp:

Mar 15, 2009, 4:56:07 PM (17 years ago)

Author:

eugene

Message:

add psastroExtract (first working version)

File:

• branches/eam_branches/eam_branch_20090312/psastro/src/psastroExtractAnalysis.c (modified) (9 diffs)

branches/eam_branches/eam_branch_20090312/psastro/src/psastroExtractAnalysis.c

@@ -20 +20 @@
 
 # define TEST_OUTPUT 1
-static pmFits *outStars = NULL;
+
+psImage *psastroExtractStar (psImage *input, float x, float y, float dX, float dY);
 
 /**
@@ -32 +33 @@
     pmReadout *readout = NULL;
     float zeropt, exptime;
+    char extname[81];
 
     // select the current recipe
@@ -40 +42 @@
     }
 
+    // force this recipe value to be false (would cause an error in psastroChooseRefstars -- no rawstars!)
+    psMetadataAddBool (recipe, PS_LIST_TAIL, "PSASTRO.MATCH.LUMFUNC", PS_META_REPLACE, "do not match luminosity functions", false);
+
     psLogMsg ("psastro", PS_LOG_INFO, "extracting bright-star subrasters");
+
+    // load the reference stars overlapping the data stars
+    psArray *refs = psastroLoadRefstars(config, "PSASTRO.EXTRACT.ASTROM");
+    if (!refs) {
+        psError (PSASTRO_ERR_UNKNOWN, false, "failed to load reference data\n");
+        return false;
+    }
+    if (refs->n == 0) {
+        psError(PSASTRO_ERR_REFSTARS, true, "no reference stars found");
+        psFree(refs);
+        return false;
+    }
+
+    if (!psastroChooseRefstars (config, refs, "PSASTRO.EXTRACT.ASTROM")) {
+        psError (PSASTRO_ERR_UNKNOWN, false, "failed to select reference data for chips\n");
+        psFree(refs);
+        return false;
+    }
+    psFree(refs);
 
     double EXTRACT_MAX_MAG = psMetadataLookupF32 (&status, recipe, "EXTRACT_MAX_MAG");
@@ -46 +70 @@
     // we have two input pmFPAfiles: PSASTRO.EXTRACT.INPUT and PSASTRO.EXTRACT.ASTROM.  both are in chip-mosaic or fpa format
     // we have already loaded the headers from PSASTRO.EXTRACT.ASTROM and determined the astrometry terms
+
+    // select the input astrometry data (also carries the refstars)
+    pmFPAfile *astrom = psMetadataLookupPtr (NULL, config->files, "PSASTRO.EXTRACT.ASTROM");
+    if (!astrom) {
+        psError(PSASTRO_ERR_CONFIG, true, "Can't find input data");
+        return false;
+    }
+    pmFPA *fpa = astrom->fpa;
 
     // select the input data sources
@@ -53 +85 @@
         return false;
     }
-    pmFPA *fpa = input->fpa;
-
-    // select the input data sources
-    pmFPAfile *astrom = psMetadataLookupPtr (NULL, config->files, "PSASTRO.EXTRACT.ASTROM");
-    if (!input) {
-        psError(PSASTRO_ERR_CONFIG, true, "Can't find input data");
-        return false;
-    }
 
     // only load data from INPUT
@@ -77 +101 @@
     EXTRACT_MAX_MAG += MagOffset;
 
+    int nExt = 0;
+
     // open the output file handle: we are just saving a series of extensions to this file
+    char *filename = psMetadataLookupStr (&status, config->arguments, "OUTPUT");
+    if (!filename) {
+      psLogMsg ("psastro", PS_LOG_INFO, "output filename not supplieda");
+      return false;
+    }
+    psFits *outStars = psFitsOpen (filename, "w");
     if (!outStars) {
-      char *filename = psMetadataLookupStr (&stats, config->arguments, "OUTPUT");
-      if (!filename) {
-        psLogMsg ("psastro", PS_LOG_INFO, "output filename not supplieda");
-        return false;
-      }
-      outStars = psFitsOpen (filename, "w");
-      if (!outStars) {
-        psError(PS_ERR_IO, false, "error opening file %s\n", filename);
-        return false;
-      }
-      // call: psFitsWriteImage (fits, NULL, subraster, 0, "extname.NN");
-    }
+      psError(PS_ERR_IO, false, "error opening file %s\n", filename);
+      return false;
+    }
+
+    pmFPAview *view = pmFPAviewAlloc (0);
 
     // open/load files as needed
@@ -109 +134 @@
           char *chipname = psMetadataLookupStr (&status, chip->concepts, "CHIP.NAME");
           psStringAppend (&filename, "bright.%s.dat", chipname);
-          FILE *f = fopen (filename, "w");
+          f = fopen (filename, "w");
           if (!f) {
             psWarning ("cannot create test output file %s\n", filename);
@@ -132 +157 @@
                 for (int i = 0; i < refstars->n; i++) {
                     pmAstromObj *ref = refstars->data[i];
-                    if (ref->Mag > REFSTAR_MASK_MAX_MAG) continue;
+                    if (ref->Mag > EXTRACT_MAX_MAG) continue;
 
                     if (TEST_OUTPUT) {
@@ -138 +163 @@
                     }
 
-                    // XXX for now, assume cell binning is 1x1 relative to chip
-                    psImage *subraster = psastroExtractStar (input->image, ref->chip->x, ref->chip->y, DX, DY);
+                    pmReadout *inReadout = pmFPAviewThisReadout (view, input->fpa);
+                    psImage *subraster = psastroExtractStar (inReadout->image, ref->chip->x, ref->chip->y, 100.0, 100.0);
+                    if (!subraster) continue;
+                    
+                    snprintf (extname, 80, "extname.%05d", nExt);
+                    psFitsWriteImage (outStars, NULL, subraster, 0, extname);
+                    nExt ++;
+                    
+                    psFree (subraster);
                 }
-                if (!pmFPAfileIOChecks (config, viewMask, PM_FPA_AFTER)) ESCAPE;
+                if (!pmFPAfileIOChecks (config, view, PM_FPA_AFTER)) ESCAPE;
             }
-            if (!pmFPAfileIOChecks (config, viewMask, PM_FPA_AFTER)) ESCAPE;
+            if (!pmFPAfileIOChecks (config, view, PM_FPA_AFTER)) ESCAPE;
         }
         if (!pmFPAfileIOChecks (config, view, PM_FPA_AFTER)) ESCAPE;
 
+        fclose (f);
     }
     if (!pmFPAfileIOChecks (config, view, PM_FPA_AFTER)) ESCAPE;
 
-    // need to close the output fits file on some basis..
-    fclose (f);
-
+    psFitsClose (outStars);
     psFree (view);
     return true;
 }
 
+psImage *psastroExtractStar (psImage *input, float x, float y, float dX, float dY) {
+
+  // skip if no pixels are on the image
+  // skip if center is not on the image
+  // if bounds fall off image, paste into a full-size image.
+
+  if (x < 0) return NULL;
+  if (y < 0) return NULL;
+  if (x >= input->numCols) return NULL;
+  if (y >= input->numRows) return NULL;
+
+  psRegion fullRegion = psRegionSet (x - dX, x + dX, y - dY, y + dY);
+  psRegion realRegion = psRegionForImage (input, fullRegion);
+  psImage *subset = psImageSubset (input, realRegion);
+
+  // if ((subset->numCols < dX) || (subset->numRows < dY)) {
+  //   psImage *fullsize = psImageAlloc (dX, dY, PS_TYPE_F32);
+  //   // fill in the fullsize image with the values from the subset
+  //   // I must already have done this elsewhere...
+  // }
+  return subset;
+}