/*****************************************************************************
******************************************************************************
**<AUTO>
**
** FILE:	stars.c
**
**<HTML>
**	This file contains C functions to select spectroscopic targets for
**	the Stars Working Group.
**</HTML>
**</AUTO>
**
**
** ENVIRONMENT:
**	ANSI C.
**
******************************************************************************
******************************************************************************
*/
#include <dervish.h>
#include "taCalibObj.h"
#include "taStars.h"
#include "arTargetType.h"
#include "arDetectionFlag.h"
#include "arObjectStatus.h"
#include "arObjectType.h"


/*****************************************************************************
******************************************************************************
**<AUTO EXTRACT>
**
** ROUTINE:	taStarsTargetSelect
**
**<HTML>
**	Set the appropriate bit in the primary target bit mask for each
**	Stars Working Group target selection criterion this object
**	satisfies.  This routine is responsible for setting the following
**	bits:
**	<menu>
**	<li> AR_TARGET_STAR_BHB
**	<li> AR_TARGET_STAR_CARBON
**	<li> AR_TARGET_STAR_BROWN_DWARF
**	<li> AR_TARGET_STAR_SUB_DWARF
**	<li> AR_TARGET_STAR_CATY_VAR
**	<li> AR_TARGET_STAR_RED_DWARF
**	<li> AR_TARGET_STAR_WHITE_DWARF
**	<li> AR_TARGET_STAR_PN
**	</menu>
**	It is assumed that the above bits are unset on entry, although other
**	bits may be set.
**</HTML>
**
** RETURNS:
**	-1	error
**	0	not a target by any criteria
**	1-255	is a target by one or more criteria --- priority returned
**		(higher number = higher priority)
**
**</AUTO>
******************************************************************************
******************************************************************************
*/
int
taStarsTargetSelect(
   TA_CALIB_OBJ *obj,               /* MOD: calibrated object to test */
   const TA_STARS_PARAMS *params    /* IN: tunable parameters*/)
{
   int priority = 0;     /* targeting priority (initialize for a non-target) */
   float u,g,r,i,z,ug,gr,ri,iz, pmDelta, hotwdug2gr;
   float dru, drg, drr, dri, drz, drug, drgr, drri, driz;
   int noblend, nointerp, accept, pmMatch;

/*
** consider only objects fainter than spectral limits in g', r', i', 
**   not too close to an edge, and with no saturated pixels
*/
   if ( obj->detection[TA_G].fiberCounts.val > params->brightLimitG &&
       obj->detection[TA_R].fiberCounts.val > params->brightLimitR &&
       obj->detection[TA_I].fiberCounts.val > params->brightLimitI &&
       !(obj->objc_flags & AR_DFLAG_SATUR) &&
       !(obj->objc_flags & AR_DFLAG_EDGE) ) {

/* exclude blended objects */
          noblend = !(obj->objc_flags & AR_DFLAG_CHILD) &&
                    !(obj->objc_flags & AR_DFLAG_BLENDED);

/* fancier code which doesn't exclude objects that photo originally
** thought were blends, then later decided they were not
**
**         noblend  = !(obj->objc_flags & AR_DFLAG_CHILD) &&
**                    ( !(obj->objc_flags & AR_DFLAG_BLENDED) ||
**                       ((obj->objc_flags & AR_DFLAG_BLENDED) &&
**                        (obj->objc_flags & AR_DFLAG_NODEBLEND)));
*/

/*
** exclude objects with interpolated pixels.
*/
            nointerp = !(obj->objc_flags & AR_DFLAG_INTERP);

            accept = noblend && nointerp;

            if (!accept) return priority;

/* use colors uncorrected for reddening for "nearby" objects */

            u = obj->detection[TA_U].psfCounts.val;
            g = obj->detection[TA_G].psfCounts.val;
            r = obj->detection[TA_R].psfCounts.val;
            i = obj->detection[TA_I].psfCounts.val;
            z = obj->detection[TA_Z].psfCounts.val;

            ug= u-g;
            gr= g-r;
            ri= r-i;
            iz= i-z;
/* 
**  use de-reddened colors (from Schlegel maps) for "halo" objects
**    or for categories where reddening would lead to stellar
**    locus contamination into selection space
*/
            dru = u - obj->reddening[TA_U];
            drg = g - obj->reddening[TA_G];
            drr = r - obj->reddening[TA_R];
            dri = i - obj->reddening[TA_I];
            drz = z - obj->reddening[TA_Z];

            drug= dru-drg;
            drgr= drg-drr;
            drri= drr-dri;
            driz= dri-drz;

            pmMatch = obj->properMotionMatch;
            pmDelta = obj->properMotionDelta;

/* Brown Dwarfs */
            if ( z > 0.0 &&
               z <= params->BWNDWRFzFaint &&
               obj->detection[TA_I].psfCounts.err <= params->BWNDWRFiErr &&
               ri >= params->BWNDWRFriBlue && iz >= params->BWNDWRFizBlue ) {
                  obj->primTarget |= AR_TARGET_STAR_BROWN_DWARF;
                  priority++;  /* Increment priority for each selection */
            }

/* Red dwarfs */
            if (obj->detection[TA_I].psfCounts.err <= params->REDDWARFiErr &&
               i <= params->REDDWARFiFaint &&
               z > 0.0 &&
               driz > params->REDDWARFizBlue && drri > params->REDDWARFriBlue){
                  obj->primTarget |= AR_TARGET_STAR_RED_DWARF;
                  priority++;
                  if (iz > 1.1) priority++;  /* more priority for reddest */
                  if (iz > 1.2) priority++;
                  if (iz > 1.4) priority++;
                  if (ri > 2.0) priority++;
                  if (ri > 2.2) priority++;
                  if (ri > 2.5) priority++;
                  if (i < 19.0) priority++;  /* and to the brightest */
            }

/* Planetary Nebulae */
            if (drgr < params->PNgrRed &&
                drri < params->PNriRed &&
                driz > params->PNizBlue &&
                drr  < params->PNrFaint &&
                drr  > params->PNrBright &&
                obj->detection[TA_G].psfCounts.err <= params->PNgErr &&
                obj->detection[TA_R].psfCounts.err <= params->PNrErr &&
                obj->detection[TA_Z].psfCounts.err <= params->PNzErr &&
                !(obj->objc_flags & AR_DFLAG_MOVED) ) {
                   obj->primTarget |= AR_TARGET_STAR_PN;
                   priority++;  /* Increment priority for each selection */
            }

/* Hot White Dwarfs */
           hotwdug2gr = drug+2*drgr;
           if(drgr < params->WHITEDWARFHOTgrRed &&
              hotwdug2gr < params->WHITEDWARFHOTuggrRed &&
              obj->detection[TA_G].psfCounts.val < params->spectralLimit &&
              params->WHITEDWARFHOTriRed > drri &&
              obj->detection[TA_G].psfCounts.val < params->spectralLimit &&
              obj->detection[TA_U].psfCounts.err <= params->WHITEDWARFHOTuErr &&
              obj->detection[TA_G].psfCounts.err <= params->WHITEDWARFHOTgErr &&
              obj->detection[TA_R].psfCounts.err <= params->WHITEDWARFHOTrErr) {
              obj->primTarget |= AR_TARGET_STAR_WHITE_DWARF;
/* set priority higher for bluer objects */
              priority += 1+(int)(.5+10.*(params->WHITEDWARFHOTuggrRed -
                                        hotwdug2gr));
              if (priority < 0) priority = 1;
              if(pmMatch == 1 && pmDelta >= params->WHITEDWARFpmLim)
                 priority++;
          }


/*
** For all other objects, require that the object_type be stellar and that
**   it has a measured r' magnitude brighter than specified limit
*/
       if ( (obj->objc_type != AR_OBJECT_TYPE_STAR) ||
            (obj->detection[TA_R].psfCounts.val <= 0.0) ||
            (obj->detection[TA_R].psfCounts.val > params->spectralLimit) )
               return priority;

/* Blue Horizontal Branch Stars */
/* funny-shaped polygon, first here's the upper rectangle */
            if ( ( (drug >= params->BHBugZ && drug < params->BHBugY) &&
                   (drgr >= params->BHBgrA && drgr < params->BHBgrB) )
/* the next lower rectangle */
              || ( (drug >= params->BHBugY && drug < params->BHBugX) &&
                   (drgr >= params->BHBgrA && drgr < params->BHBgrC) )
/* the top part of the irregular pentagon */
              || ( (drug >= params->BHBugX && drug < params->BHBugW) &&
                   (drgr >= params->BHBgrA &&
                    drgr <= params->BHBgrD+(params->BHBgrE-params->BHBgrD)*
                       (drug-params->BHBugX)/(params->BHBugW-params->BHBugX)))
/* finally, the bottom part of the irregular pentagon */
              || ( (drug >= params->BHBugW && drug <= params->BHBugV) &&
                   (drgr >= params->BHBgrA &&
                    drgr <= params->BHBgrE-(params->BHBgrE-params->BHBgrD)*
                       (drug-params->BHBugW)/(params->BHBugV-params->BHBugW)))){
                obj->primTarget |= AR_TARGET_STAR_BHB;
                priority++;  /* Increment priority for each selection */
                if (drgr < 0.1)  priority++;  /* preference for bluer stars in box */
                if (drgr < 0.0)  priority++;  /* more preference for even bluer stars in box */
                if (drg  < 19.0) priority++; /* preference for brighter stars */
             }

/* Carbon Stars */
             if ( (drgr >= params->CARBONgrBlue && drri >= params->CARBONriBlue
                && driz >= params->CARBONizBlue
                && (drri <= params->CARBONriIntercept +
                    params->CARBONriSlope*drgr) )
                || (drgr >= params->CARBONgrDusty) )  {
                  obj->primTarget |= AR_TARGET_STAR_CARBON;
                  priority++;  /* Increment priority for each selection */
                  if(drgr >= params->CARBONgrRed)
                     priority++;
         }

/* Low-luminosity subdwarfs */
            if(drgr > params->LOWLUMSUBDWRFgrBlue &&
               drri > params->LOWLUMSUBDWRFriBlue &&
               drri < params->LOWLUMSUBDWRFriRed  &&
            obj->detection[TA_G].psfCounts.err <= params->LOWLUMSUBDWRFgErr &&
            obj->detection[TA_R].psfCounts.err <= params->LOWLUMSUBDWRFrErr &&
            obj->detection[TA_I].psfCounts.err <= params->LOWLUMSUBDWRFiErr) {
                  obj->primTarget |= AR_TARGET_STAR_SUB_DWARF;
                  priority++;  /* Increment priority for each selection */
                  if(pmMatch == 1 && pmDelta >= params->LOWLUMSUBDWRFpmLim)
                     priority++; /* Increment priority if it has sig. p.m. */
            }

/* Cataclysmic variables */
            if(g <= params->CATYVARgFaint &&
               ug <= params->CATYVARugRed && 
               gr <= params->CATYVARgrRed &&
               ri >= params->CATYVARriBlue &&
               iz >= params->CATYVARizBlue ) {
                  obj->primTarget |= AR_TARGET_STAR_CATY_VAR;
                  priority++;  /* Increment priority for each selection */
            }

/* Low luminosity white dwarfs **
**    there are two subclasses of white dwarfs, hot (above) and low-luminosity
**    the low-luminosity white dwarfs always get higher priority than do
**    the hot white dwarfs
*/
           if( (drug > params->WHITEDWARFugBlue ||
                dru  > params->WHITEDWARFuBright) &&
                obj->detection[TA_G].psfCounts.err <= params->WHITEDWARFgErr &&
                obj->detection[TA_R].psfCounts.err <= params->WHITEDWARFrErr &&
                obj->detection[TA_I].psfCounts.err <= params->WHITEDWARFiErr &&
                obj->detection[TA_Z].psfCounts.err <= params->WHITEDWARFzErr &&
                drgr > params->WHITEDWARFgrBlue &&
                drgr < params->WHITEDWARFgrRed &&
                driz < params->WHITEDWARFizRed &&
                drgr > params->WHITEDWARFgrSlope*drri +
                   params->WHITEDWARFgrIntercept &&
                params->WHITEDWARFrizBlue < params->WHITEDWARFriFactor*drri+
                params->WHITEDWARFizFactor*driz  &&
                params->WHITEDWARFrizRed  > params->WHITEDWARFriFactor*drri+
                params->WHITEDWARFizFactor*driz ) {
                   obj->primTarget |= AR_TARGET_STAR_WHITE_DWARF;
                   priority+=127;   /* this is how we set low-lum wd priorities
                                  ** to be higher than hot wd */
                   if(pmMatch == 1 && pmDelta >= params->WHITEDWARFpmLim)
                      priority++;
           }
   }
/* Successful return */
   return priority;
}
