Added CONCENTRATION_MODEL and finalized CLASS_ST parameters.

Version number pushed to 2.8.4.

configure

 #! /bin/sh
 # Guess values for system-dependent variables and create Makefiles.
-# Generated by GNU Autoconf 2.61 for sextractor 2.8.3.
+# Generated by GNU Autoconf 2.61 for sextractor 2.8.4.
 #
 # Report bugs to <bertin@iap.fr>.
 #
@@ -728,8 +728,8 @@ SHELL=${CONFIG_SHELL-/bin/sh}
 # Identity of this package.
 PACKAGE_NAME='sextractor'
 PACKAGE_TARNAME='sextractor'
-PACKAGE_VERSION='2.8.3'
-PACKAGE_STRING='sextractor 2.8.3'
+PACKAGE_VERSION='2.8.4'
+PACKAGE_STRING='sextractor 2.8.4'
 PACKAGE_BUGREPORT='bertin@iap.fr'
 ac_unique_file="src/makeit.c"
@@ -1405,7 +1405,7 @@ if test "$ac_init_help" = "long"; then
   # Omit some internal or obsolete options to make the list less imposing.
   # This message is too long to be a string in the A/UX 3.1 sh.
   cat <<_ACEOF
-\`configure' configures sextractor 2.8.3 to adapt to many kinds of systems.
+\`configure' configures sextractor 2.8.4 to adapt to many kinds of systems.
 Usage: $0 [OPTION]... [VAR=VALUE]...
@@ -1475,7 +1475,7 @@ fi
 if test -n "$ac_init_help"; then
   case $ac_init_help in
-     short | recursive ) echo "Configuration of sextractor 2.8.3:";;
+     short | recursive ) echo "Configuration of sextractor 2.8.4:";;
    esac
   cat <<\_ACEOF
@@ -1599,7 +1599,7 @@ fi
 test -n "$ac_init_help" && exit $ac_status
 if $ac_init_version; then
   cat <<\_ACEOF
-sextractor configure 2.8.3
+sextractor configure 2.8.4
 generated by GNU Autoconf 2.61
 Copyright (C) 1992, 1993, 1994, 1995, 1996, 1998, 1999, 2000, 2001,
@@ -1613,7 +1613,7 @@ cat >config.log <<_ACEOF
 This file contains any messages produced by compilers while
 running configure, to aid debugging if configure makes a mistake.
-It was created by sextractor $as_me 2.8.3, which was
+It was created by sextractor $as_me 2.8.4, which was
 generated by GNU Autoconf 2.61.  Invocation command line was
   $ $0 $@
@@ -2290,7 +2290,7 @@ fi
 # Define the identity of the package.
  PACKAGE='sextractor'
- VERSION='2.8.3'
+ VERSION='2.8.4'
 cat >>confdefs.h <<_ACEOF
@@ -27297,7 +27297,7 @@ exec 6>&1
 # report actual input values of CONFIG_FILES etc. instead of their
 # values after options handling.
 ac_log="
-This file was extended by sextractor $as_me 2.8.3, which was
+This file was extended by sextractor $as_me 2.8.4, which was
 generated by GNU Autoconf 2.61.  Invocation command line was
   CONFIG_FILES    = $CONFIG_FILES
@@ -27350,7 +27350,7 @@ Report bugs to <bug-autoconf@gnu.org>."
 _ACEOF
 cat >>$CONFIG_STATUS <<_ACEOF
 ac_cs_version="\\
-sextractor config.status 2.8.3
+sextractor config.status 2.8.4
 configured by $0, generated by GNU Autoconf 2.61,
   with options \\"`echo "$ac_configure_args" | sed 's/^ //; s/[\\""\`\$]/\\\\&/g'`\\"

configure.ac

@@ -6,7 +6,7 @@ define([AC_CACHE_LOAD],)
 define([AC_CACHE_SAVE],)
 
 # This is your standard Bertin source code...
-AC_INIT(sextractor, 2.8.3, [bertin@iap.fr])
+AC_INIT(sextractor, 2.8.4, [bertin@iap.fr])
 AC_CONFIG_SRCDIR(src/makeit.c)
 AC_CONFIG_AUX_DIR(autoconf)
 AM_CONFIG_HEADER(config.h)

man/sex.1

-.TH SEXTRACTOR "1" "February 2009" "SWarp 2.8.3" "User Commands"
+.TH SEXTRACTOR "1" "March 2009" "SWarp 2.8.4" "User Commands"
 .SH NAME
 sex \- extract a source catalog from an astronomical FITS image
 .SH SYNOPSIS

src/catout.c

@@ -9,7 +9,7 @@
 *
 *	Contents:	functions for output of catalog data.
 *
-*	Last modify:	25/09/2008
+*	Last modify:	18/03/2009
 *
 *%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 */
@@ -205,7 +205,9 @@ void	updateparamflags()
   FLAG(obj2.x_prof) |= FLAG(obj2.y_prof)
 			| FLAG(obj2.xw_prof)
 			| FLAG(obj2.poserra_prof)
-			| FLAG(obj2.poserrcxx_prof);
+			| FLAG(obj2.poserrcxx_prof)
+			| FLAG(obj2.prof_concentration)
+			| FLAG(obj2.prof_class_star);
   FLAG(obj2.mag_prof) |= FLAG(obj2.magerr_prof);
   FLAG(obj2.flux_prof) |= FLAG(obj2.mag_prof) | FLAG(obj2.fluxerr_prof);
   FLAG(obj2.prof_spheroid_mag) |= FLAG(obj2.prof_spheroid_magerr);

src/paramprofit.h

@@ -9,7 +9,7 @@
 *
 *	Contents:	Model-fitting parameter list for catalog data.
 *
-*	Last modify:	25/09/2008
+*	Last modify:	18/03/2009
 *
 *%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 */
@@ -167,6 +167,13 @@
 	&outobj2.prof_eps2, H_FLOAT, T_FLOAT, "%10.6f", "",
 	"src.ellipticity;stat.fit;instr.det", ""},
 
+  {"CONCENTRATION_MODEL", "Concentration parameter from model-fitting",
+	&outobj2.prof_concentration, H_FLOAT, T_FLOAT, "%7.4f", "",
+	"src.morph.param", ""},
+  {"CLASS_STAR_MODEL", "S/G classifier from model-fitting",
+	&outobj2.prof_class_star, H_FLOAT, T_FLOAT, "%7.4f", "",
+	"src.class.starGalaxy", ""},
+
   {"FLUX_BACKOFFSET", "Background offset from fitting",
 	&outobj2.prof_offset_flux, H_FLOAT, T_FLOAT, "%12.7g", "count",
 	"instr.skyLevel;arith.diff;stat.fit.param", "ct"},

src/pattern.c

@@ -9,7 +9,7 @@
 *
 *	Contents:	Generate and handle image patterns for image fitting.
 *
-*	Last modify:	20/11/2008
+*	Last modify:	18/03/2009
 *
 *%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 */
@@ -140,7 +140,7 @@ INPUT	Pointer to pattern structure.
 OUTPUT	-.
 NOTES	-.
 AUTHOR	E. Bertin (IAP)
-VERSION	01/10/2008
+VERSION	18/03/2009
  ***/
 void	pattern_fit(patternstruct *pattern, profitstruct *profit)
   {
@@ -162,7 +162,7 @@ void	pattern_fit(patternstruct *pattern, profitstruct *profit)
   for (p=0; p<nvec; p++)
     {
     profit_convolve(profit, inpix);
-    profit_resample(profit, inpix, outpix);
+    profit_resample(profit, inpix, outpix, 1.0);
     outpix1 = pattern->lmodpix;
     for (p2=0; p2<=p; p2++)
       {

src/profit.c

@@ -9,7 +9,7 @@
 *
 *	Contents:	Fit an arbitrary profile combination to a detection.
 *
-*	Last modify:	29/01/2009
+*	Last modify:	18/03/2009
 *
 *%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 */
@@ -168,17 +168,18 @@ OUTPUT	Pointer to an allocated fit structure (containing details about the
 	fit).
 NOTES	It is a modified version of the lm_minimize() of lmfit.
 AUTHOR	E. Bertin (IAP)
-VERSION	29/01/2009
+VERSION	18/03/2009
  ***/
 void	profit_fit(profitstruct *profit,
 		picstruct *field, picstruct *wfield,
 		objstruct *obj, obj2struct *obj2)
   {
+    profitstruct	pprofit;
     patternstruct *pattern;
     psfstruct		*psf;
     checkstruct		*check;
     double		*oldparaminit,
-			psf_fwhm, oldchi2, a , cp,sp, emx2,emy2,emxy;
+			psf_fwhm, oldchi2, a , cp,sp, emx2,emy2,emxy, dchi2;
     int			i,j,p, oldniter, nparam, ncomp;
 
   nparam = profit->nparam;
@@ -251,10 +252,6 @@ void	profit_fit(profitstruct *profit,
 
   profit_resetparams(profit);
 
-//for (p=0; p<nparam; p++)
-//printf("%g ", profit->paraminit[p]);
-//printf("\n");
-
 the_gal++;
 
 /* Actual minimisation */
@@ -287,11 +284,6 @@ the_gal++;
   for (p=0; p<profit->nparam; p++)
     profit->param[p] = profit->paraminit[p];
 
-//printf("--> ");
-//for (p=0; p<profit->nparam; p++)
-//printf("%g     %g %g\n", profit->param[p], profit->parammin[p], profit->parammax[p]);
-//printf("(%d)\n", profit->niter);
-
 /* CHECK-Images */
   if ((check = prefs.check[CHECK_SUBPROFILES]))
     {
@@ -305,7 +297,6 @@ the_gal++;
     addcheck(check, profit->lmodpix, profit->objnaxisn[0],profit->objnaxisn[1],
 		profit->ix,profit->iy, 1.0);
     }
-
 /* Fill measurement parameters */
   if (FLAG(obj2.prof_vector))
     {
@@ -326,11 +317,24 @@ the_gal++;
     i = profit->paramindex[PARAM_X];
     j = profit->paramindex[PARAM_Y];
 /*-- Model coordinates follow the FITS convention (first pixel at 1,1) */
+    if (profit->paramlist[PARAM_X])
+      {
       obj2->x_prof = profit->ix + *profit->paramlist[PARAM_X] + 1.0;
-    obj2->y_prof = profit->iy + *profit->paramlist[PARAM_Y] + 1.0;
       obj2->poserrmx2_prof = emx2 = profit->covar[i*(nparam+1)];
+      }
+    else
+      emx2 = 0.0;
+    if (profit->paramlist[PARAM_Y])
+      {
+      obj2->y_prof = profit->iy + *profit->paramlist[PARAM_Y] + 1.0;
       obj2->poserrmy2_prof = emy2 = profit->covar[j*(nparam+1)];
+      }
+    else
+      emy2 = 0.0;
+    if (profit->paramlist[PARAM_X] && profit->paramlist[PARAM_Y])
       obj2->poserrmxy_prof = emxy = profit->covar[i+j*nparam];
+    else
+      emxy = 0.0;
 
 /*-- Error ellipse parameters */
     if (FLAG(obj2.poserra_prof))
@@ -508,6 +512,38 @@ the_gal++;
       }
     }
 
+  if (FLAG(obj2.prof_class_star))
+    {
+    pprofit = *profit;
+    memset(pprofit.paramindex, 0, PARAM_NPARAM*sizeof(int));
+    memset(pprofit.paramlist, 0, PARAM_NPARAM*sizeof(double *));
+    pprofit.nparam = 0;
+    QMALLOC(pprofit.prof, profstruct *, 1);
+    pprofit.prof[0] = prof_init(&pprofit, PROF_DIRAC);
+    QMALLOC(pprofit.covar, double, pprofit.nparam*pprofit.nparam);
+    pprofit.nprof = 1;
+    profit_resetparams(&pprofit);
+    if (profit->paramlist[PARAM_X] && profit->paramlist[PARAM_Y])
+      {
+      pprofit.paraminit[pprofit.paramindex[PARAM_X]] = *profit->paramlist[PARAM_X];
+      pprofit.paraminit[pprofit.paramindex[PARAM_Y]] = *profit->paramlist[PARAM_Y];
+      }
+    pprofit.paraminit[pprofit.paramindex[PARAM_DISK_FLUX]] = profit->flux;
+    pprofit.niter = profit_minimize(&pprofit, PROFIT_MAXITER);
+    dchi2 = 0.5*(pprofit.chi2 - obj2->prof_chi2);
+    obj2->prof_class_star = dchi2 < 50.0?
+	(dchi2 > -50.0? 2.0/(1.0+exp(dchi2)) : 2.0) : 0.0;
+    if (profit->flux > 0.0 && pprofit.flux > 0.0)
+      obj2->prof_concentration = -2.5*log10(pprofit.flux / profit->flux);
+    else  if (profit->flux > 0.0)
+      obj2->prof_concentration = 99.0;
+    else  if (pprofit.flux > 0.0)
+      obj2->prof_concentration = -99.0;
+    prof_end(pprofit.prof[0]);
+    free(pprofit.prof);
+    free(pprofit.covar);
+    }
+
 /* clean up. */
   free(profit->modpix);
   free(profit->psfpix);
@@ -517,9 +553,7 @@ the_gal++;
   free(profit->objweight);
   free(profit->resi);
   free(oldparaminit);
-//for (ix=0; ix<profit->nparam;ix++)
-//printf("%10.5f ", sqrt(profit->covar[ix*profit->nparam+ix]));
-//printf("\n");
+
   return;
   }
 
@@ -739,7 +773,7 @@ INPUT	Profile-fitting structure,
 OUTPUT	Vector of residuals.
 NOTES	-.
 AUTHOR	E. Bertin (IAP)
-VERSION	15/09/2008
+VERSION	18/03/2009
  ***/
 double	*profit_residuals(profitstruct *profit, picstruct *field,
 		picstruct *wfield, double *param, double *resi)
@@ -750,10 +784,16 @@ double	*profit_residuals(profitstruct *profit, picstruct *field,
 	profit->modnaxisn[0]*profit->modnaxisn[1]*sizeof(double));
   for (p=0; p<profit->nparam; p++)
     profit->param[p] = param[p];
+/* Simple PSF shortcut */
+  if (profit->nprof == 1 && profit->prof[0]->code == PROF_DIRAC)
+    profit_resample(profit, profit->psfpix, profit->lmodpix, *profit->prof[0]->flux);
+  else
+    {
     for (p=0; p<profit->nprof; p++)
       prof_add(profit->prof[p], profit);
     profit_convolve(profit, profit->modpix);
-  profit_resample(profit, profit->modpix, profit->lmodpix);
+    profit_resample(profit, profit->modpix, profit->lmodpix, 1.0);
+    }
   profit_compresi(profit, resi);
 
   return resi;
@@ -770,7 +810,7 @@ INPUT	Profile-fitting structure,
 OUTPUT	Vector of residuals.
 NOTES	-.
 AUTHOR	E. Bertin (IAP)
-VERSION	18/09/2008
+VERSION	11/03/2009
  ***/
 double	*profit_compresi(profitstruct *profit, double *resi)
   {
@@ -781,7 +821,6 @@ double	*profit_compresi(profitstruct *profit, double *resi)
    int		npix, ix1,ix2;
   
 /* Compute vector of residuals */
-  npix = profit->objnaxisn[0]*profit->objnaxisn[1];
   resit = resi;
   objpix = profit->objpix;
   objweight = profit->objweight;
@@ -791,6 +830,7 @@ double	*profit_compresi(profitstruct *profit, double *resi)
   x1c = (double)(profit->objnaxisn[0]/2);
   rmin = profit->obj2->hl_radius / 2.0;
   x2 = -(double)(profit->objnaxisn[1]/2);
+  npix = -profit->nparam;
   for (ix2=profit->objnaxisn[1]; ix2--; x2+=1.0)
     {
     x1 = -x1c;
@@ -804,11 +844,14 @@ double	*profit_compresi(profitstruct *profit, double *resi)
         *(resit++) = val2*PROFIT_DYNPARAM;
 //        *(resit++) = val2*(rmin/(sqrt(r2)+rmin));
         error += val2*val2;
+        npix++;
         }
       }
     }
 
-  profit->chi2 = PROFIT_DYNPARAM*PROFIT_DYNPARAM*error/npix;
+  if (npix <= 0)
+    npix = 1;
+  profit->chi2 = PROFIT_DYNPARAM*PROFIT_DYNPARAM*error;
 
   return resi;
   }
@@ -822,9 +865,10 @@ INPUT	Profile-fitting structure.
 OUTPUT	-.
 NOTES	-.
 AUTHOR	E. Bertin (IAP)
-VERSION	15/09/2008
+VERSION	18/03/2009
  ***/
-void	profit_resample(profitstruct *profit, double *inpix, PIXTYPE *outpix)
+void	profit_resample(profitstruct *profit, double *inpix, PIXTYPE *outpix,
+	double factor)
   {
    double	posin[2], posout[2], dnaxisn[2],
 		*dx,*dy,
@@ -854,7 +898,7 @@ void	profit_resample(profitstruct *profit, double *inpix, PIXTYPE *outpix)
     {
     posin[0] = (posout[0] - xcout)*invpixstep + xcin;
     posin[1] = (posout[1] - ycout)*invpixstep + ycin;
-    flux += ((*(outpix++) = (PIXTYPE)(interpolate_pix(posin, inpix,
+    flux += ((*(outpix++) = (PIXTYPE)(factor*interpolate_pix(posin, inpix,
 		profit->modnaxisn, INTERP_LANCZOS3))));
     for (d=0; d<2; d++)
       if ((posout[d]+=1.0) < dnaxisn[d])
@@ -1569,7 +1613,7 @@ INPUT	Pointer to the profit structure,
 	Upper boundary to the parameter.
 OUTPUT	RETURN_OK if the parameter is registered, RETURN_ERROR otherwise.
 AUTHOR	E. Bertin (IAP)
-VERSION	16/04/2008
+VERSION	15/03/2009
  ***/
 int	profit_setparam(profitstruct *profit, paramenum paramtype,
 		double param, double parammin, double parammax)
@@ -1580,7 +1624,7 @@ int	profit_setparam(profitstruct *profit, paramenum paramtype,
 /* Check whether the parameter has already be registered */
   if ((paramptr=profit->paramlist[(int)paramtype]))
     {
-    index = paramptr - profit->param;
+    index = profit->paramindex[(int)paramtype];
     profit->paraminit[index] = param;
     profit->parammin[index] = parammin;
     profit->parammax[index] = parammax;
@@ -1804,6 +1848,14 @@ profstruct	*prof_init(profitstruct *profit, proftypenum profcode)
       profit_addparam(profit, PARAM_OUTRING_FLUX, &prof->flux);
       profit_addparam(profit, PARAM_OUTRING_WIDTH, &prof->featwidth);
       break;
+    case PROF_DIRAC:
+      prof->naxis = 2;
+      prof->pix = NULL;
+      prof->typscale = 1.0;
+      profit_addparam(profit, PARAM_X, &prof->x[0]);
+      profit_addparam(profit, PARAM_Y, &prof->x[1]);
+      profit_addparam(profit, PARAM_DISK_FLUX, &prof->flux);
+      break;
     case PROF_SERSIC_TABEX:	/* An example of tabulated profile */
       prof->naxis = 3;
       width =  prof->naxisn[0] = PROFIT_PROFRES;
@@ -1928,7 +1980,9 @@ void	prof_add(profstruct *prof, profitstruct *profit)
 
   scaling = profit->pixstep / prof->typscale;
 
-/* Compute Profile CD matrix */
+  if (prof->code!=PROF_DIRAC)
+    {
+/*-- Compute Profile CD matrix */
     ctheta = cos(*prof->posangle*DEG);
     stheta = sin(*prof->posangle*DEG);
     saspect = fabs(*prof->aspect);
@@ -1940,6 +1994,7 @@ void	prof_add(profstruct *prof, profitstruct *profit)
     cd12 = xscale*stheta;
     cd21 = -yscale*stheta;
     cd22 = yscale*ctheta;
+    }
 
   dx1 = 0.0;	/* Shifting operations have been moved to profit_resample() */
   dx2 = 0.0;	/* Shifting operations have been moved to profit_resample() */
@@ -2203,6 +2258,12 @@ width = 3.0;
           }
         }
       break;
+    case PROF_DIRAC:
+      memset(profit->pmodpix, 0,
+	profit->modnaxisn[0]*profit->modnaxisn[1]*sizeof(float));
+      profit->pmodpix[profit->modnaxisn[0]/2
+		+ (profit->modnaxisn[1]/2)*profit->modnaxisn[0]] = 1.0;
+      break;
     default:
 /*---- Tabulated profile: remap each pixel */
 /*---- Initialize multi-dimensional counters */

src/profit.h

@@ -9,7 +9,7 @@
 *
 *	Contents:	Include file for profit.c.
 *
-*	Last modify:	03/10/2008
+*	Last modify:	18/03/2009
 *
 *%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 */
@@ -46,7 +46,7 @@ One must have:	PROFIT_NITER > 0
 
 typedef enum		{PROF_BACK, PROF_SERSIC, PROF_DEVAUCOULEURS,
 			PROF_EXPONENTIAL, PROF_ARMS, PROF_BAR, PROF_INRING,
-			PROF_OUTRING, PROF_SERSIC_TABEX, PROF_NPROF}
+			PROF_OUTRING, PROF_SERSIC_TABEX, PROF_DIRAC, PROF_NPROF}
 				proftypenum; /* Profile code */
 
 typedef enum	{INTERP_NEARESTNEIGHBOUR, INTERP_BILINEAR, INTERP_LANCZOS2,
@@ -172,7 +172,7 @@ void		prof_add(profstruct *prof, profitstruct *profit),
 			void *data, int iflag, int iter, int nfev ),
 		profit_psf(profitstruct *profit),
 		profit_resample(profitstruct *profit, double *inpix,
-			PIXTYPE *outpix),
+			PIXTYPE *outpix, double factor),
 		profit_resetparam(profitstruct *profit, paramenum paramtype),
 		profit_resetparams(profitstruct *profit),
 		profit_unboundtobound(profitstruct *profit, double *param);

src/types.h

@@ -9,7 +9,7 @@
 *
 *	Contents:	global type definitions.
 *
-*	Last modify:	25/09/2008
+*	Last modify:	18/03/2009
 *
 *%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 */
@@ -334,6 +334,8 @@ typedef struct
   double	prof_mx2w, prof_my2w, prof_mxyw;/* WORLD profile model moments*/
   float		prof_eps1, prof_eps2;		/* Profile model ellip.vector */
   float		prof_e1, prof_e2;		/* Profile model ellip.vector */
+  float		prof_class_star;		/* Model-fitting star/gal class*/
+  float		prof_concentration;		/* Model-fitting concentration*/
   float		prof_offset_flux;		/* Background offset */
   float		prof_offset_fluxerr;		/* RMS error */
   float		prof_spheroid_flux;		/* Spheroid total flux */