/* @(#)fit_orders.c	17.1.1.1 (ESO-IPG) 01/25/02 17:51:58 */
/* 

   Otmar Stahl

   fit_orders.c 

   fit order positions with polynomial, each order individually

*/

/* system includes */

#include <stdio.h>

/* general Midas includes */

#include <midas_def.h>

/* FEROS specific includes */

#include <proto_nrutil.h>
#include <proto_mutil.h>
#include <echelle.h>

#ifdef __STDC__
int fit_orders 
(
 int tid, int nact, int fit_deg, float tab[], int icol[],
 float rnull
 )
#else
     int fit_orders 
(
 tid, nact, fit_deg, tab, icol, rnull
 )
     int tid, nact, fit_deg, icol[];
     float tab[], rnull;
#endif

{
  int ncol, nrow, nsort, allcol, allrow, dunit, i, j, jj, null[7];
  double xpos[MAXLINES], ypos[MAXLINES], xposs, yposs;
  double *ad;
  char line[80];

  int first_usable_order, last_usable_order;
  ad = dvector (1, fit_deg+1);

  first_usable_order = 0;
  last_usable_order = 0;
  TCIGET(tid, &ncol, &nrow, &nsort, &allcol, &allrow);
  SCDWRI(tid, "FITORD", &fit_deg, 1, 1, &dunit);
  
  for (j = 1; j <= nact; j++)
    {
      jj = 0;
      for (i = 1; i <= nrow; i++)
	{
	  TCRRDR(tid, i, 5, icol, tab, null);
	  /* select entries for order from table */
	  if(((int) tab[2] == j) && tab[0] != rnull)
	    {
	      ypos[jj] = (double) tab[0];
	      xpos[jj] = (double) tab[1];
	      jj++;
	    }
	}

      /* fit positions */

      sprintf(line, "order %i of %i", j, nact);
      SCTPUT(line);
      if (jj >= fit_deg + 1)
	{
	  if (first_usable_order == 0) first_usable_order = j;
	  last_usable_order = j;
	  lfit(xpos, ypos, jj, ad, fit_deg, dpoly);

	  /* write descriptors */

	  sprintf (line,"FIT%04i", j);
	  SCDWRD(tid, line, &ad[1], 1, fit_deg, &dunit);
	  for (i = 1; i <= nrow; i++)
	    {
	      TCRRDR (tid, i, 5, icol, tab, null);
	      if(((int) tab[2] == j) && tab[0] != rnull)

		/* write fitted positions and residuals to output table */

		{
		  xposs = (double) tab[1];
		  yposs = (double) tab[0];
		  tab[4] = (float) eval_dpoly (xposs, ad, fit_deg);
		  tab[5] = yposs-tab[4];
		  TCEWRR(tid, i, icol[4], &tab[4]);
		  TCEWRR(tid, i, icol[5], &tab[5]);
		}
	    }
	}
      else
	{
	  if (last_usable_order < first_usable_order)
	    last_usable_order = first_usable_order;

	  /* must be set; otherwise there will be problems */

	  sprintf(line, "order %i cannot be fitted", j);
	  SCTPUT(line);

	  /* no fit can be done */

	}
    }
  SCDWRI(tid, "FIRSTORD", &first_usable_order, 1, 1, &dunit);
  SCDWRI(tid, "ECHORD", &last_usable_order, 1, 1, &dunit);
  free_dvector (ad, 1, fit_deg+1);
  return 0;
}