Index: trunk/psastro/src/psastroLoadGhosts.c
===================================================================
--- trunk/psastro/src/psastroLoadGhosts.c	(revision 41367)
+++ trunk/psastro/src/psastroLoadGhosts.c	(revision 41434)
@@ -41,4 +41,10 @@
     }
 
+// To select good stars, I want them to have a PSFmodel fit (PM_SOURCE_MODE_PSFMODEL & PM_SOURCE_MODE_FITTED) and exlcude things that are CR, EXT, FAIL, POOR 
+// Stars are allowed to be saturated, as long as the PSF model fits
+# define PHOT_SOURCE_MASK (PM_SOURCE_MODE_FAIL | PM_SOURCE_MODE_BLEND | PM_SOURCE_MODE_BADPSF | \
+                           PM_SOURCE_MODE_DEFECT | PM_SOURCE_MODE_CR_LIMIT | PM_SOURCE_MODE_EXT_LIMIT | \
+                           PM_SOURCE_MODE_POOR) // Mask to apply to sources for rejection
+
 /**
  * calculate ghost FPA and Chip positions for the stars loaded on the FPA
@@ -101,5 +107,5 @@
     pmCell *cell = NULL;
     pmReadout *readout = NULL;
-    float zeropt, exptime, MAX_MAG;
+    float zeropt, exptime,MAX_MAG, INSTR_MAX_MAG;
     psMetadata *md = NULL;
     psPolynomial2D *centerX = NULL;
@@ -139,5 +145,5 @@
 
     // raise an error if the config is broken
-    if (!psastroZeroPointFromRecipe (&zeropt, &exptime, &MAX_MAG, NULL, fpa, recipe)) {
+    if (!psastroZeroPointFromRecipe (&zeropt, &exptime, &INSTR_MAX_MAG, NULL, fpa, recipe)) {
         psError(PSASTRO_ERR_CONFIG, true, "failed to load zeropt data from recipe");
 	goto escape;
@@ -147,5 +153,5 @@
     // use the zero point and exposure time to convert to apparent mags: M_ap = M_inst + C_0 + 2.5*log(exptime)
     float MagOffset = zeropt + 2.5*log10(exptime);
-    MAX_MAG += MagOffset;
+    MAX_MAG = INSTR_MAX_MAG + MagOffset;
 
     // this loop selects the matched stars for all chips
@@ -163,4 +169,119 @@
                 if (! readout->data_exists) { continue; }
 
+                //First, select bright sources for ghosts from the detections. Later on, do the same using the refcat
+                // That way we can include things like movers creating crosstalk features.
+                psArray *calstars = psMetadataLookupPtr (&status, readout->analysis, "PSASTRO.CALSTARS");
+                if (calstars == NULL) { continue; }
+
+                // identify the bright stars of interest
+                for (int i = 0; i < calstars->n; i++) {
+                    pmAstromObj *cal = calstars->data[i];
+		    psTrace("psastro.ghost",5,"Begin ghost %d/%ld: MAX_MAG: %g; ref_mag: %g @ (%.10g,%.10g)",i,calstars->n,INSTR_MAX_MAG,cal->Mag,cal->sky->r * 180 / M_PI,cal->sky->d * 180.0 / M_PI);
+
+                    if (cal->Mag > INSTR_MAX_MAG) {
+                        continue;
+                    }
+
+		    psastroGhost *ghost = psastroGhostAlloc ();
+		    ghost->srcFP->x = cal->FP->x; 
+		    ghost->srcFP->y = cal->FP->y;
+
+		    double rSrc = hypot (cal->FP->x, cal->FP->y);
+     	            double theta0 = atan2(cal->FP->x,cal->FP->y);
+
+                    if(mirCheck) {
+                         //TdB: first mirror the reference star positions (around the 0,0 pixel) using the radial offset coefficients and the ghost/star angle
+		        double ghost_offset_rad = psPolynomial1DEval (mirrorRad, rSrc);
+		        double ghost_x_fpa_mirror = cal->FP->x + ((cal->FP->x*-1.)/abs(cal->FP->x)*abs(cos(theta0)*ghost_offset_rad));
+		        double ghost_y_fpa_mirror = cal->FP->y + ((cal->FP->y*-1.)/abs(cal->FP->y)*abs(sin(theta0)*ghost_offset_rad));
+
+		        // Now use the mirrored position together with the 2D ghost center fitting to get the actual ghost position in FPA coords 
+		        ghost->FP->x = ghost_x_fpa_mirror + psPolynomial2DEval(centerX, ghost_x_fpa_mirror, ghost_y_fpa_mirror);
+		        ghost->FP->y = ghost_y_fpa_mirror + psPolynomial2DEval(centerY, ghost_x_fpa_mirror, ghost_y_fpa_mirror);
+                    }  
+                    if(!mirCheck) {
+                        //Use the old-style ghost position determination
+                        ghost->FP->x = -cal->FP->x + psPolynomial2DEval(centerX, -cal->FP->x, -cal->FP->y);
+                        ghost->FP->y = -cal->FP->y + psPolynomial2DEval(centerY, -cal->FP->x, -cal->FP->y);
+                    }
+
+		    ghost->inner.major = psPolynomial1DEval (innerMajor, rSrc);
+		    ghost->inner.minor = psPolynomial1DEval (innerMinor, rSrc);
+		    ghost->inner.theta = atan2(cal->FP->y, cal->FP->x);
+
+		    ghost->outer.major = psPolynomial1DEval (outerMajor, rSrc);
+		    ghost->outer.minor = psPolynomial1DEval (outerMinor, rSrc);
+		    ghost->outer.theta = atan2(cal->FP->y, cal->FP->x);
+
+		    // report instrumental ghost star mags
+		    ghost->Mag = cal->Mag;
+		    
+		    // XXX this code yields a single chip: we need to provide results for any chips
+		    // which encompass the full size of the ghost
+		    pmChip *ghostChip = psastroFindChip (&ghost->chip->x, &ghost->chip->y, fpa, -ghost->srcFP->x, -ghost->srcFP->y);
+/* 		    float star_chip_x = ghost->chip->x; */
+/* 		    float star_chip_y = ghost->chip->y; */
+		    ghostChip = psastroFindChip (&ghost->chip->x, &ghost->chip->y, fpa, ghost->FP->x, ghost->FP->y);
+                    //psLogMsg ("psastro", PS_LOG_INFO, "-> DETEC model chip position: %f, %f %f, %f %g\n", ghost->chip->x, ghost->chip->y,ghost->srcFP->x, ghost->srcFP->y,cal->Mag+MagOffset);
+
+		    if (ghostChip) {
+		      psTrace("psastro.ghost",5,"   in ghost: %d/%ld: ref: (%f,%f) or (%s) ghost: (%f,%f) or (%f,%f,%s) %f inner: (%f,%f,%f) outer: (%f,%f,%f)",
+			      i,calstars->n,
+			      cal->FP->x,cal->FP->y,
+			      psMetadataLookupStr(NULL,chip->concepts,"CHIP.NAME"),
+			      ghost->FP->x,ghost->FP->y,
+			      ghost->chip->x,ghost->chip->y,psMetadataLookupStr(NULL,ghostChip->concepts,"CHIP.NAME"),
+			      rSrc,
+			      ghost->inner.major,ghost->inner.minor,ghost->inner.theta,
+			      ghost->outer.major,ghost->outer.minor,ghost->outer.theta);
+/* 		      psWarning("   GHOST_DATA: %d/%ld: ref: (%f,%f) or (%f,%f,%s) ghost: (%f,%f) or (%f,%f,%s) %f inner: (%f,%f,%f) outer: (%f,%f,%f)", */
+/* 			      i,calstars->n, */
+/* 			      cal->FP->x,cal->FP->y, */
+/* 			      ref_chip_x,ref_chip_y,psMetadataLookupStr(NULL,chip->concepts,"CHIP.NAME"), */
+/* 			      ghost->FP->x,ghost->FP->y, */
+/* 			      ghost->chip->x,ghost->chip->y,psMetadataLookupStr(NULL,ghostChip->concepts,"CHIP.NAME"), */
+/* 			      rSrc, */
+/* 			      ghost->inner.major,ghost->inner.minor,ghost->inner.theta, */
+/* 			      ghost->outer.major,ghost->outer.minor,ghost->outer.theta); */
+		    }
+		    else {
+		      psTrace("psastro.ghost",5,"   in ghost: %d/%ld: ref: (%f,%f) ghost: (%f,%f) or (%f,%f,%s) inner: (%f,%f,%f) outer: (%f,%f,%f)",
+			      i,calstars->n,
+			      cal->FP->x,cal->FP->y,
+			      ghost->FP->x,ghost->FP->y,
+			      ghost->chip->x,ghost->chip->y,"NONE",
+			      ghost->inner.major,ghost->inner.minor,ghost->inner.theta,
+			      ghost->outer.major,ghost->outer.minor,ghost->outer.theta);
+		    }		      
+		      			    
+		    if (!ghostChip) goto skip;
+		    if (!ghostChip->cells) goto skip;
+		    if (!ghostChip->cells->n) goto skip;
+
+		    
+		    pmCell *ghostCell = ghostChip->cells->data[0];
+		    if (!ghostCell) goto skip;
+		    if (!ghostCell->readouts) goto skip;
+		    if (!ghostCell->readouts->n) goto skip;
+		    pmReadout *ghostReadout = ghostCell->readouts->data[0];
+		    if (!ghostReadout) goto skip;
+
+		    psArray *ghosts = psMetadataLookupPtr (&status, ghostReadout->analysis, "PSASTRO.GHOSTS");
+		    if (ghosts == NULL) { 
+			ghosts = psArrayAllocEmpty (100);
+			if (!psMetadataAdd (ghostReadout->analysis, PS_LIST_TAIL, "PSASTRO.GHOSTS", PS_DATA_ARRAY, "astrometry matches", ghosts)) {
+			  psError(PSASTRO_ERR_CONFIG, false, "failure to add ghosts to readout");
+			  goto escape;
+			}
+			psFree (ghosts);
+		    }
+
+		    psArrayAdd (ghosts, 100, ghost);
+
+		skip:
+		    
+		    psFree (ghost);
+                }
+
                 // select the raw objects for this readout (loaded in psastroChooseRefstars.c)
 		// XXX : note that we place limits on the refstar sample in psastroChooseRefstars.c:
@@ -173,4 +294,5 @@
                     pmAstromObj *ref = refstars->data[i];
 		    psTrace("psastro.ghost",5,"Begin ghost %d/%ld: MAX_MAG: %g; ref_mag: %g @ (%.10g,%.10g)",i,refstars->n,MAX_MAG,ref->Mag,ref->sky->r * 180 / M_PI,ref->sky->d * 180.0 / M_PI);
+
                     if (ref->Mag > MAX_MAG) continue;
 
@@ -224,6 +346,5 @@
 /* 		    float ref_chip_y = ghost->chip->y; */
 		    ghostChip = psastroFindChip (&ghost->chip->x, &ghost->chip->y, fpa, ghost->FP->x, ghost->FP->y);
-		    // fprintf (stderr, "-> model chip position: %f, %f\n", ghost->chip->x, ghost->chip->y);
-
+                    //psLogMsg ("psastro", PS_LOG_INFO, "-> model chip position: %f, %f %f, %f %g\n", ghost->chip->x, ghost->chip->y,ghost->srcFP->x, ghost->srcFP->y,ref->Mag);
 
 		    if (ghostChip) {
@@ -257,15 +378,15 @@
 		    }		      
 		      			    
-		    if (!ghostChip) goto skip;
-		    if (!ghostChip->cells) goto skip;
-		    if (!ghostChip->cells->n) goto skip;
+		    if (!ghostChip) goto skipref;
+		    if (!ghostChip->cells) goto skipref;
+		    if (!ghostChip->cells->n) goto skipref;
 
 		    
 		    pmCell *ghostCell = ghostChip->cells->data[0];
-		    if (!ghostCell) goto skip;
-		    if (!ghostCell->readouts) goto skip;
-		    if (!ghostCell->readouts->n) goto skip;
+		    if (!ghostCell) goto skipref;
+		    if (!ghostCell->readouts) goto skipref;
+		    if (!ghostCell->readouts->n) goto skipref;
 		    pmReadout *ghostReadout = ghostCell->readouts->data[0];
-		    if (!ghostReadout) goto skip;
+		    if (!ghostReadout) goto skipref;
 
 		    psArray *ghosts = psMetadataLookupPtr (&status, ghostReadout->analysis, "PSASTRO.GHOSTS");
@@ -281,8 +402,9 @@
 		    psArrayAdd (ghosts, 100, ghost);
 
-		skip:
+		skipref:
 		    
 		    psFree (ghost);
                 }
+
             }
         }