/*      FILE: /mxtools/src/psfstar/img_star.c
 *   PURPOSE: 
 *    AUTHOR: Kenneth J. Mighell (mighell@noao.edu)
 *  LANGUAGE: ANSI C
 *      DATE: 2001NOV07
 * COPYRIGHT: (C) 2001 Assoc. of Universities for Research in Astronomy Inc.
 */

#include <mx.h>
#include "inc.h"
#include <assert.h>

int inc_ImageS_Add_Star_i7(
  struct mxip_image_s *image,
  double se_p,
  double sx_p,
  double sy_p,
  struct mxip_image_s *psf,
  int PSF_CPXJ_p,
  int PSF_CPYJ_p
) {

  char mxfunc[] = "inc_ImageS_Add_Star_i7";
  int mxstatus = 0;
  int status = 1;
  int debug = MX_TRUE;
#define DEBUG (debug==MX_FALSE)

  double se;
  int sxj;
  int syj;

  int inx;
  int iny;
  int pnx;
  int pny;
  int snx;
  int sny;

  double tpxd;
  double ipxd;
  double ipxd0;
  double ppxd;

  int sllxj;
  int sllyj;
  int surxj;
  int suryj;

  int illxj;
  int illyj;
  int iurxj;
  int iuryj;

  int deltaxj;
  int deltayj;

  int ixj;
  int iyj;
  int pxj;
  int pyj;

  inx = image->nxi;
  iny = image->nyi;

  pnx = psf->nxi;
  pny = psf->nyi;

  se = se_p;
  sxj = sx_p - 1;
  syj = sy_p - 1;

  illxj = image->llx0i;
  illyj = image->lly0i;
  iurxj = illxj + (image->nxi) - 1;
  iuryj = illyj + (image->nyi) - 1;

  deltaxj = sxj - PSF_CPXJ_p;
  deltayj = syj - PSF_CPYJ_p;
  sllxj = deltaxj - (psf->llx0i);
  sllyj = deltayj - (psf->lly0i);
  surxj = sllxj + (psf->nxi) - 1;
  suryj = sllyj + (psf->nyi) - 1;

  /* Is STAR LLX right of IMAGE URX ? */
  if (sllxj>iurxj) {
    if (DEBUG) printf("sllxj>iurxj!\n");
    goto ok;
  }

  /* Is STAR URX left of IMAGE LLX? */
  if (surxj<illxj) {
    if (DEBUG) printf("surxj<illxj!\n");
    goto ok;
  }

  /* Is STAR LLY above IMAGE URY ? */
  if (sllyj>iuryj) {
    if (DEBUG) printf("sllyj>iuryj!\n");
    goto ok;
  }

  /* Is STAR URY below IMAGE LLY ? */
  if (suryj<illyj) {
    if (DEBUG) printf("suryj<illyj!\n");
    goto ok;
  }

/*
 * At least one pixel needs to be changed...
 */

  /* Is STAR LLX left of IMAGE LLX ? */
  if (sllxj<illxj) {
    sllxj = illxj;
    if (DEBUG) printf("sllxj<illxj: sllxj is now %d\n", sllxj);
  }

  /* Is STAR URX right of IMAGE URX ? */
  if (surxj>iurxj) {
    surxj = iurxj;
    if (DEBUG) printf("surxj>iurxj: surxj is now %d\n", surxj);
  }

  /* Is STAR LLY below IMAGE LLY ? */
  if (sllyj<illyj) {
    sllyj = illyj;
    if (DEBUG) printf("sllyj<illyj: sllyj is now %d\n", sllyj);
  }

  /* Is STAR URY above IMAGE URY ? */
  if (suryj>iuryj) {
    suryj = iuryj;
    if (DEBUG) printf("suryj>iuryj: suryj is now %d\n", suryj);
  }

  /* Number of pixels in the X direction */
  snx = surxj - sllxj + 1;
  if (DEBUG) printf("snx is %d\n",snx);

  /* Do we have at least one pixel in the X direction? */
  if (snx<1) {
    if (DEBUG) printf("snx<1!\n");
    goto ok;
  }

  /* Do we have no more than PNX pixels in the X dirction? */
  if (snx>pnx) {
    if (DEBUG) printf("snx>pnx!\n");
    goto ok;
  }

  /* Number of pixels in the Y direction */
  sny = suryj - sllyj + 1;
  if (DEBUG) printf("sny is %d\n",sny);

  /* Do we have at least one pixel in the Y direction? */
  if (sny>pny) {
    if (DEBUG) printf("sny>pny!\n");
    goto ok;
  }

  /* Do we have no more than PNY pixels in the Y dirction? */
  if (sny>pny) {
    if (DEBUG) printf("sny>pny!\n");
    goto ok;
  }

  if (DEBUG) printf(" :-) \n");

  for (iyj=sllyj; iyj<=suryj; ++iyj) {
    assert((iyj>=0)&&(iyj<iny));   /* What, me worry? */
    pyj = iyj - deltayj;
    assert((pyj>=0)&&(pyj<pny));   /* What, me worry? */
    for (ixj=sllxj; ixj<=surxj; ++ixj) {
      assert((ixj>=0)&&(ixj<inx)); /* What, me worry? */
      pxj = ixj - deltaxj;
      assert((pxj>=0)&&(pxj<pnx)); /* What, me worry? */
      ipxd0 = image->matrixd[iyj][ixj];
      ppxd = psf->matrixd[pyj][pxj];
      if (ppxd!=0.0) {
	tpxd = ppxd * se; /* stellar PSF * stellar intensity */
        ipxd = ipxd0 + tpxd;
        image->matrixd[iyj][ixj] = ipxd; /* add stellar pixel to image */
	if (DEBUG) printf(
	  "%6d %6d %10g -> %10g | %6d %6d %10g -> %10g\n", 
          pyj, pxj, ppxd, tpxd, iyj, ixj, ipxd0, ipxd
        );
      }
    }
  }

ok:
  status = 0;
  goto bye;

mx_error1:
  mxp_tmpmsg_init_f0();
mx_error2:
  mxp_errmsg_append_f3 (mxfunc, mxstatus, MX.tmpmsg);
  goto bye;

bye:
  return(status);

}
#undef DEBUG
/* end-of-file */