Changeset 24691 for trunk/magic/remove/src/streaksremove.c

Timestamp:

Jul 6, 2009, 4:24:20 PM (17 years ago)

Author:

bills

Message:

when destreaking warp and diff skyfiles, censor any sources whose coordinates are masked

File:

• trunk/magic/remove/src/streaksremove.c (modified) (6 diffs)

trunk/magic/remove/src/streaksremove.c

@@ -17 +17 @@
 static void writeImages(streakFiles *sf, bool exciseImageCube);
 static void updateAstrometry(streakFiles *sfiles);
+static void censorSources(streakFiles *sfiles, psU32 maskStreak);
 
 int
@@ -188 +189 @@
         }
 
+        censorSources(sfiles, maskStreak);
+
         // write the destreaked "temporary" images and the recovery images
         writeImages(sfiles, exciseImageCube);
+
 
         psLogMsg("streaksremove", PS_LOG_INFO, "time to process component %d: %f\n", sfiles->extnum, psTimerClear("PROCESS_COMPONENT"));
@@ -259 +263 @@
     fprintf(stderr, "\t-weight WEIGHT.fits: weight file to de-streak\n");
     fprintf(stderr, "\t-replace: replace the input images with the output\n");
-    fprintf(stderr, "\t-remove: remove the original image after processing (requires -replace)\n");
     fprintf(stderr, "\t-keepnonwarped: do not exise pixels that were not part of difference processing\n");
     fprintf(stderr, "\t-transparent val: instead of setting excicsed pixel to NAN add val\n");
@@ -335 +338 @@
     bool gotReplace = false;
     if ((argnum = psArgumentGet(argc, argv, "-replace"))) {
+        if (!nebulousImage) {
+            psError(PS_ERR_UNKNOWN, true, "replace is only supported for nebulous files");
+            usage();
+        }
         gotReplace = true;
         psArgumentRemove(argnum, &argc, argv);
@@ -416 +423 @@
         psMetadataAddStr(config->arguments, PS_LIST_TAIL, "INPUT.WEIGHT", 0,
                 "name of input weight image", argv[argnum]);
+        psArgumentRemove(argnum, &argc, argv);
+    }
+    if ((argnum = psArgumentGet(argc, argv, "-sources"))) {
+        psArgumentRemove(argnum, &argc, argv);
+        bool nebulousSources = IN_NEBULOUS(argv[argnum]);
+        if (nebulousSources != nebulousImage) {
+            psError(PS_ERR_UNKNOWN, true, "sources file must have %snebulous path with %s image path\n",
+                nebulousImage ? "" : "non ", nebulousImage ? "nebulous" : "regular");
+            usage();
+        }
+        psMetadataAddStr(config->arguments, PS_LIST_TAIL, "SOURCES", 0,
+                "name of input sources file", argv[argnum]);
         psArgumentRemove(argnum, &argc, argv);
     }
@@ -757 +776 @@
 }
 
+// read a sources file (.cmf) and remove any rows whose coordinate is convered by a
+// streak mask
+static void
+censorSources(streakFiles *sfiles, psU32 maskStreak)
+{
+    if ((!sfiles->inSources) || (!sfiles->outMask)) {
+        return;
+    }
+    psImage *maskImage = sfiles->outMask->image;
+    if (!maskImage) {
+        psError(PS_ERR_IO, false, "maskImage is null!");
+        streaksExit("", PS_EXIT_PROG_ERROR);
+    }
+
+    sFile *in = sfiles->inSources;
+    sFile *out = sfiles->outSources;
+
+    in->header = psFitsReadHeader(NULL, in->fits);
+    if (!in->header) {
+        psError(PS_ERR_IO, false, "failed to read header from %s", in->resolved_name);
+        streaksExit("", PS_EXIT_DATA_ERROR);
+    }
+
+    bool status;
+    psString extname = psMetadataLookupStr(&status, in->header, "EXTNAME");
+    if (!extname) {
+        psError(PS_ERR_IO, false, "failed to find extname in header of %s", in->resolved_name);
+        streaksExit("", PS_EXIT_DATA_ERROR);
+    }
+
+    psArray *inTable = psFitsReadTable(in->fits);
+    if (!inTable->n) {
+        psError(PS_ERR_IO, false, "table in %s is empty", in->resolved_name);
+        streaksExit("", PS_EXIT_DATA_ERROR);
+    }
+
+    psArray *outTable = psArrayAllocEmpty(inTable->n);
+    int j = 0;
+    int numCensored = 0;
+    for (int i = 0 ; i < inTable->n; i++) {
+        psMetadata *row = inTable->data[i];
+
+        psF32 x = psMetadataLookupF32 (&status, row, "X_PSF");
+        psF32 y = psMetadataLookupF32 (&status, row, "Y_PSF");
+        
+        psU32 mask = psImageGet(maskImage, x, y);
+
+        // Key the source if the center pixel is not masked with maskStreak
+        if (! (mask & maskStreak) ) {
+            psArraySet(outTable, j++, row);
+        } else {
+            numCensored++;
+        }
+    }
+
+    printf("Censored %d sources\n", numCensored);
+
+    // get rid of unused elements (don't know if this is necessary)
+    psArrayRealloc(outTable, j);
+
+    addDestreakKeyword(in->header);
+    if (! psFitsWriteTable(out->fits, in->header, outTable, extname)) {
+        psError(PS_ERR_IO, false, "failed to write table to %s", out->resolved_name);
+        streaksExit("", PS_EXIT_DATA_ERROR);
+    }
+
+    if (!psFitsClose(out->fits)) {
+        psError(PS_ERR_IO, false, "failed to close table %s", out->resolved_name);
+        streaksExit("", PS_EXIT_DATA_ERROR);
+    }
+}