Changeset 19200 for trunk/pswarp/src/pswarpTransformSources.c

Timestamp:

Aug 25, 2008, 3:04:38 PM (18 years ago)

Author:

Paul Price

Message:

Propagate more source characteristics, esp the mode (specifies CR/EXT).

File:

• trunk/pswarp/src/pswarpTransformSources.c (modified) (2 diffs)

trunk/pswarp/src/pswarpTransformSources.c

@@ -24 +24 @@
     psArray *outSources = psMemIncrRefCounter(psMetadataLookupPtr(&mdok, output->analysis, "PSPHOT.SOURCES")); // Target sources
     if (!outSources) {
-        outSources = psArrayAllocEmpty(SOURCE_ARRAY_BUFFER);
-        psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "PSPHOT.SOURCES", PS_DATA_ARRAY,
-                         "Warped sources", outSources);
+        outSources = psArrayAllocEmpty(SOURCE_ARRAY_BUFFER);
+        psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "PSPHOT.SOURCES", PS_DATA_ARRAY,
+                         "Warped sources", outSources);
     }
 
@@ -33 +33 @@
     // this does not cost us so much time.
     for (int i = 0; i < inSources->n; i++) {
-        pmSource *source = inSources->data[i]; // Source of interest
-        pmModel *model = source->modelPSF; // Model for this source
-        float xIn, yIn;             // Coordinates of source
-        xIn = model->params->data.F32[PM_PAR_XPOS] - input->image->col0;
-        yIn = model->params->data.F32[PM_PAR_YPOS] - input->image->row0;
+        pmSource *source = inSources->data[i]; // Source of interest
+        pmModel *model = source->modelPSF; // Model for this source
+        float xIn, yIn;             // Coordinates of source
+        xIn = model->params->data.F32[PM_PAR_XPOS] - input->image->col0;
+        yIn = model->params->data.F32[PM_PAR_YPOS] - input->image->row0;
 
-        int xGrid, yGrid;           // Grid coordinates for local map
-        if (!pswarpMapGridSetGrid(sourceGrid, xIn + 0.5, yIn + 0.5, &xGrid, &yGrid)) {
-            psError(PS_ERR_UNKNOWN, false, "Unable to get grid coordinates for source at %f,%f\n",
-                    xIn, yIn);
-            psFree(outSources);
-            psFree(sourceGrid);
-            return false;
-        }
-        if (xGrid < 0 || xGrid >= sourceGrid->nXpts || yGrid < 0 || yGrid >= sourceGrid->nYpts) {
-            // It's not even on the grid
-            // XXX how can this happen?
-            continue;
-        }
+        int xGrid, yGrid;           // Grid coordinates for local map
+        if (!pswarpMapGridSetGrid(sourceGrid, xIn + 0.5, yIn + 0.5, &xGrid, &yGrid)) {
+            psError(PS_ERR_UNKNOWN, false, "Unable to get grid coordinates for source at %f,%f\n",
+                    xIn, yIn);
+            psFree(outSources);
+            psFree(sourceGrid);
+            return false;
+        }
+        if (xGrid < 0 || xGrid >= sourceGrid->nXpts || yGrid < 0 || yGrid >= sourceGrid->nYpts) {
+            // It's not even on the grid
+            // XXX how can this happen?
+            continue;
+        }
 
-        pswarpMap *map = sourceGrid->maps[xGrid][yGrid]; // Locally linear transformation
-        double xOut, yOut;          // Output coordinates
-        if (!pswarpMapApply(&xOut, &yOut, map, xIn + 0.5, yIn + 0.5)) {
-            psError(PS_ERR_UNKNOWN, false, "Unable to transform coordinates for source at %f,%f\n",
-                    xIn, yIn);
-            psFree(outSources);
-            psFree(sourceGrid);
-            return false;
-        }
-        xOut += output->image->col0 - 0.5;
-        yOut += output->image->row0 - 0.5;
-        if (xOut < minX || xOut > maxX || yOut < minY || yOut > maxY) {
-            // It's not in the output image
-            continue;
-        }
+        pswarpMap *map = sourceGrid->maps[xGrid][yGrid]; // Locally linear transformation
+        double xOut, yOut;          // Output coordinates
+        if (!pswarpMapApply(&xOut, &yOut, map, xIn + 0.5, yIn + 0.5)) {
+            psError(PS_ERR_UNKNOWN, false, "Unable to transform coordinates for source at %f,%f\n",
+                    xIn, yIn);
+            psFree(outSources);
+            psFree(sourceGrid);
+            return false;
+        }
+        xOut += output->image->col0 - 0.5;
+        yOut += output->image->row0 - 0.5;
+        if (xOut < minX || xOut > maxX || yOut < minY || yOut > maxY) {
+            // It's not in the output image
+            continue;
+        }
 
-        // Generate the new source in the output frame
-        pmSource *new = pmSourceAlloc(); // New source
-        new->peak = pmPeakAlloc(xOut, yOut, source->peak->flux, PM_PEAK_LONE);
-        new->peak->flux = source->peak->flux;
-        new->type = PM_SOURCE_TYPE_STAR;
-        new->psfMag = source->psfMag;
-        new->errMag = source->errMag;
-        new->sky = source->sky;
-        new->skyErr = source->skyErr;
-        new->pixWeight = source->pixWeight;
-        new->modelPSF = pmModelAlloc(source->modelPSF->type);
+        // Generate the new source in the output frame
+        //
+        // Magnitudes will be off if there's any change in scale, but for our purposes (mainly x,y and
+        // relative flux) that's OK.
+        pmSource *new = pmSourceAlloc(); // New source
+        new->peak = pmPeakAlloc(xOut, yOut, source->peak->flux, PM_PEAK_LONE);
+        new->peak->flux = source->peak->flux;
+        new->type = source->type;
+        new->mode = source->mode;
+        new->psfMag = source->psfMag;
+        new->extMag = source->extMag;
+        new->errMag = source->errMag;
+        new->apMag = source->apMag;
+        new->pixWeight = source->pixWeight;
+        new->psfChisq = source->psfChisq;
+        new->crNsigma = source->crNsigma;
+        new->extNsigma = source->extNsigma;
+        new->sky = source->sky;
+        new->skyErr = source->skyErr;
+
+        new->modelPSF = pmModelAlloc(source->modelPSF->type);
 
 #if 0
-        // XXX Note that this will not set the correct axes
-        pmPSF_AxesToModel(new->modelPSF->params->data.F32,
-                          pmPSF_ModelToAxes(source->modelPSF->params->data.F32, 20.0));
+        // XXX Note that this will not set the correct axes
+        pmPSF_AxesToModel(new->modelPSF->params->data.F32,
+                          pmPSF_ModelToAxes(source->modelPSF->params->data.F32, 20.0));
 #endif
 
-        // Propagate the position erorrs
-        float dxIn, dyIn;           // Errors in input coordinates
-        dxIn = model->dparams->data.F32[PM_PAR_XPOS];
-        dyIn = model->dparams->data.F32[PM_PAR_YPOS];
+        // Propagate the position erorrs
+        float dxIn, dyIn;           // Errors in input coordinates
+        dxIn = model->dparams->data.F32[PM_PAR_XPOS];
+        dyIn = model->dparams->data.F32[PM_PAR_YPOS];
 
-        float dxOut, dyOut;         // Errors in output coordinates
-        dxOut = sqrt(PS_SQR(map->Xx * dxIn) + PS_SQR(map->Xy * dyIn));
-        dyOut = sqrt(PS_SQR(map->Yx * dxIn) + PS_SQR(map->Yy * dyIn));
+        float dxOut, dyOut;         // Errors in output coordinates
+        dxOut = sqrt(PS_SQR(map->Xx * dxIn) + PS_SQR(map->Xy * dyIn));
+        dyOut = sqrt(PS_SQR(map->Yx * dxIn) + PS_SQR(map->Yy * dyIn));
 
-        new->modelPSF->params->data.F32[PM_PAR_XPOS] = xOut;
-        new->modelPSF->params->data.F32[PM_PAR_YPOS] = yOut;
-        new->modelPSF->dparams->data.F32[PM_PAR_XPOS] = dxOut;
-        new->modelPSF->dparams->data.F32[PM_PAR_YPOS] = dyOut;
+        new->modelPSF->params->data.F32[PM_PAR_XPOS] = xOut;
+        new->modelPSF->params->data.F32[PM_PAR_YPOS] = yOut;
+        new->modelPSF->dparams->data.F32[PM_PAR_XPOS] = dxOut;
+        new->modelPSF->dparams->data.F32[PM_PAR_YPOS] = dyOut;
 
-        psArrayAdd(outSources, SOURCE_ARRAY_BUFFER, new);
-        psFree(new);                // Drop reference
+        psArrayAdd(outSources, SOURCE_ARRAY_BUFFER, new);
+        psFree(new);                // Drop reference
     }
     psFree(sourceGrid);