Changeset 32397

Timestamp:

Sep 11, 2011, 4:31:22 PM (15 years ago)

Author:

bills

Message:

ppMops memory fix

Location:

tags/ipp-20110622/ppTranslate/src

Files:

• . (modified) (1 prop)
• ppMopsRead.c (modified) (2 diffs)

tags/ipp-20110622/ppTranslate/src/ppMopsRead.c

@@ -95 +95 @@
         psFree(header);
 
+#ifndef USE_OLD_READ_TABLE
+        // use new fits table functions
+
+        psFitsTable *table = psFitsReadTableNew(fits); // Table of interest
+        if (!table) {
+            psError(PS_ERR_IO, false, "Unable to read table %d", i);
+            return false;
+        }
+        psFitsClose(fits);
+
+        double plateScale = 0.0;        // Plate scale
+        long numGood = 0;               // Number of good rows
+        for (long row = 0; row < size; row++) {
+
+            psU32 flags = psFitsTableGetU32(NULL, table, row, "FLAGS");
+            if (flags & SOURCE_MASK) {
+                psTrace("ppMops.read", 10, "Discarding row %ld from input %d because of flags: %ud", row, i, flags);
+                continue;
+            }
+
+            det->x->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "X_PSF");
+            det->y->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "Y_PSF");
+            det->ra->data.F64[numGood] = DEG_TO_RAD(psFitsTableGetF64(NULL, table, row, "RA_PSF"));
+            det->dec->data.F64[numGood] = DEG_TO_RAD(psFitsTableGetF64(NULL, table, row, "DEC_PSF"));
+            det->mag->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_INST_MAG");
+            det->magErr->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_INST_MAG_SIG");
+            det->chi2->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_CHISQ");
+            det->dof->data.S32[numGood] = psFitsTableGetS32(NULL, table, row, "PSF_NDOF");
+            det->cr->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "CR_NSIGMA");
+            det->extended->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "EXT_NSIGMA");
+            det->psfMajor->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_MAJOR");
+            det->psfMinor->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_MINOR");
+            det->psfTheta->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_THETA");
+            det->quality->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_QF");
+            det->numPix->data.S32[numGood] = psFitsTableGetS32(NULL, table, row, "PSF_NPIX");
+            det->xxMoment->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "MOMENTS_XX");
+            det->xyMoment->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "MOMENTS_XY");
+            det->yyMoment->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "MOMENTS_YY");
+            det->flags->data.U32[numGood] = psFitsTableGetU32(NULL, table, row, "FLAGS");
+            det->diffSkyfileId->data.S64[numGood] = diffSkyfileId;
+            det->naxis1->data.S32[numGood] = naxis1;
+            det->naxis2->data.S32[numGood] = naxis2;
+            det->nPos->data.S32[numGood] = psFitsTableGetS32(NULL, table, row, "DIFF_NPOS");
+            det->fPos->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "DIFF_FRATIO");
+            det->nRatioBad->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "DIFF_NRATIO_BAD");
+            det->nRatioMask->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "DIFF_NRATIO_MASK");
+            det->nRatioAll->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "DIFF_NRATIO_ALL");
+
+            //Additions of 2010-10-25
+            if (args->version == 2) {
+              //Values are set only if the version is 2
+              if (skyChipPsfVersion == 2) {
+                det->psfInstFlux->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_INST_FLUX");
+                det->psfInstFluxSig->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_INST_FLUX_SIG");
+                det->apMag->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "AP_MAG");
+                det->apMagRaw->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "AP_MAG_RAW");
+                det->apMagRadius->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "AP_MAG_RADIUS");
+                det->apFlux->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "AP_FLUX");
+                det->apFluxSig->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "AP_FLUX_SIG");
+                det->peakFluxAsMag->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PEAK_FLUX_AS_MAG");
+                det->calPsfMag->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "CAL_PSF_MAG");
+                det->calPsfMagSig->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "CAL_PSF_MAG_SIG");
+                det->sky->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "SKY");
+                det->skySig->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "SKY_SIGMA");
+                det->qualityPerfect->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "PSF_QF_PERFECT");
+                det->momentsR1->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "MOMENTS_R1");
+                det->momentsRH->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "MOMENTS_RH");
+                det->kronFlux->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "KRON_FLUX");
+                det->kronFluxErr->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "KRON_FLUX_ERR");
+                det->kronFluxInner->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "KRON_FLUX_INNER");
+                det->kronFluxOuter->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "KRON_FLUX_OUTER");
+                det->diffRP->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "DIFF_R_P");
+                det->diffSnP->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "DIFF_SN_P");
+                det->diffRM->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "DIFF_R_M");
+                det->diffSnM->data.F32[numGood] = psFitsTableGetF32(NULL, table, row, "DIFF_SN_M");
+                det->flags2->data.U32[numGood] = psFitsTableGetU32(NULL, table, row, "FLAGS2");
+                det->ippIdet->data.U32[numGood] = psFitsTableGetU32(NULL, table, row, "IPP_IDET");
+                det->nFrames->data.U16[numGood] = psFitsTableGetU16(NULL, table, row, "N_FRAMES");
+                det->padding->data.S16[numGood] = psFitsTableGetS16(NULL, table, row, "PADDING");
+              } else {
+                det->psfInstFlux->data.F32[numGood] = NAN;
+                det->psfInstFluxSig->data.F32[numGood] = NAN;
+                det->apMag->data.F32[numGood] = NAN;
+                det->apMagRaw->data.F32[numGood] = NAN;
+                det->apMagRadius->data.F32[numGood] = NAN;
+                det->apFlux->data.F32[numGood] = NAN;
+                det->apFluxSig->data.F32[numGood] = NAN;
+                det->peakFluxAsMag->data.F32[numGood] = NAN;
+                det->calPsfMag->data.F32[numGood] = NAN;
+                det->calPsfMagSig->data.F32[numGood] = NAN;
+                det->sky->data.F32[numGood] = NAN;
+                det->skySig->data.F32[numGood] = NAN;
+                det->qualityPerfect->data.F32[numGood] = NAN;
+                det->momentsR1->data.F32[numGood] = NAN;
+                det->momentsRH->data.F32[numGood] = NAN;
+                det->kronFlux->data.F32[numGood] = NAN;
+                det->kronFluxErr->data.F32[numGood] = NAN;
+                det->kronFluxInner->data.F32[numGood] = NAN;
+                det->kronFluxOuter->data.F32[numGood] = NAN;
+                det->diffRP->data.F32[numGood] = NAN;
+                det->diffSnP->data.F32[numGood] = NAN;
+                det->diffRM->data.F32[numGood] = NAN;
+                det->diffSnM->data.F32[numGood] = NAN;
+                det->flags2->data.U32[numGood] = 0;
+                det->ippIdet->data.U32[numGood] = 0;
+                det->nFrames->data.U16[numGood] = 0;
+                det->padding->data.S16[numGood] = 0;
+              }
+            }
+
+            // Calculate error in RA, Dec
+            double xErr = psFitsTableGetF64(NULL, table, row, "X_PSF_SIG"); //SC: Warning! Promotion of F32
+            double yErr = psFitsTableGetF64(NULL, table, row, "Y_PSF_SIG"); //SC: Warning! Promotion of F32
+            double scale = psFitsTableGetF64(NULL, table, row, "PLTSCALE"); //SC: Warning! Promotion of F32
+            double angle = psFitsTableGetF64(NULL, table, row, "POSANGLE"); //SC: Warning! Promotion of F32
+
+            if (!isfinite(det->x->data.F32[numGood]) || !isfinite(det->y->data.F32[numGood]) ||
+                !isfinite(det->ra->data.F64[numGood]) || !isfinite(det->dec->data.F64[numGood]) ||
+                !isfinite(det->mag->data.F32[numGood]) || !isfinite(det->magErr->data.F32[numGood]) ||
+                !isfinite(xErr) || !isfinite(yErr) || !isfinite(scale) || !isfinite(angle)) {
+                psTrace("ppMops.read", 10,
+                        "Discarding row %ld from input %d because of non-finite values: "
+                        "%f %f %lf %lf %f %f %f %f %f %f",
+                        row, i,
+                        det->x->data.F32[numGood], det->y->data.F32[numGood],
+                        det->ra->data.F64[numGood], det->dec->data.F64[numGood],
+                        det->mag->data.F32[numGood], det->magErr->data.F32[numGood],
+                        xErr, yErr, scale, angle);
+                continue;
+            }
+
+            // XXX Not at all sure I've got the angles around the right way here...
+            double cosAngle = cos(angle), sinAngle = sin(angle);
+            double cosAngle2 = PS_SQR(cosAngle), sinAngle2 = PS_SQR(sinAngle);
+            double xErr2 = PS_SQR(xErr), yErr2 = PS_SQR(yErr);
+            double errScale = scale / 3600.0;
+            det->raErr->data.F64[numGood] = errScale * sqrt(cosAngle2 * xErr2 + sinAngle2 * yErr2);
+            det->decErr->data.F64[numGood] = errScale * sqrt(sinAngle2 * xErr2 + cosAngle2 * yErr2);
+
+            det->mask->data.U8[numGood] = 0;
+            plateScale += scale;
+            numGood++;
+        }
+#else
+    // OLD way of reading fits tables read it into an array of metadata objects, one for each row. Uses lots and lots
+    // and lots of memory.
         psArray *table = psFitsReadTable(fits); // Table of interest
         if (!table) {
@@ -236 +382 @@
             numGood++;
         }
+#endif
         det->seeing *= ((float) plateScale) / ((float) numGood);