/* @(#)fitsinf.c	17.1.1.1 (ESO-DMD) 01/25/02 17:39:09 */
/*===========================================================================
  Copyright (C) 1995 European Southern Observatory (ESO)
 
  This program is free software; you can redistribute it and/or 
  modify it under the terms of the GNU General Public License as 
  published by the Free Software Foundation; either version 2 of 
  the License, or (at your option) any later version.
 
  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.
 
  You should have received a copy of the GNU General Public 
  License along with this program; if not, write to the Free 
  Software Foundation, Inc., 675 Massachusetts Ave, Cambridge, 
  MA 02139, USA.
 
  Correspondence concerning ESO-MIDAS should be addressed as follows:
	Internet e-mail: midas@eso.org
	Postal address: European Southern Observatory
			Data Management Division 
			Karl-Schwarzschild-Strasse 2
			D 85748 Garching bei Muenchen 
			GERMANY
===========================================================================*/

/*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
.COPYRIGHT  (c)  1995   European Southern Observatory
.IDENT      fitsinf.c
.LAUGUAGE   C
.AUTHOR     P.Grosbol   ESO/IPG
.KEYWORDS   MIDAS frames, information
.COMMENT    get informations from MIDAS frames
.VERSION    1.0  1991-Mar-23 : Creation,   PJG
.VERSION    1.1  1991-May-15 : Check for SCFGET return,   PJG
.VERSION    1.2  1991-Jun-03 : Always set bfac and bzero,   PJG
.VERSION    1.3  1991-Sep-24 : Get NULL's for tables,   PJG
.VERSION    1.4  1992-Aug-07 : Correct format for C*n table columns, PJG
.VERSION    1.5  1993-Apr-12 : Correct MIN/MAX computation, PJG
.VERSION    1.6  1993-Sep-17 : Correct fitstbl because of TCBGET
.VERSION    1.7  1993-Oct-26 : Update to new TC + prototypes, PJG
.VERSION    1.8  1995-Jan-23 : Do not update LHCUTS, correct check, PJG
.VERSION    1.9  1996-Jun-24 : Exchange MAXLONG with 2147483645.0, PJG
.VERSION    2.0  1998-Mar-09 : use existing BSCALE, BZERO if FITS file already
.VERSION    2.1  1998-Apr-22 : Use I2 format of I1 table columns, PJG
---------------------------------------------------------------------*/
#include   <osparms.h>
#include   <computer.h>
#include   <fitsdef.h>
#include   <fitskwt.h>
#include   <fileexts.h>
/*
#include   <midas_def.h>
*/
#include   <tblsys.h>

#define    BSIZE       10240         /* size of internal buffer      */
#define    MXLB           81         /* max. char. in line buffer    */

static     int         fpexc;        /* fp- exception flag           */

static     SDEF         sdef;        /* Image scaling information    */
static     TXDEF       txdef;        /* table definitions            */
static     FDEF    fdef[MXF];        /* information on columns       */

SDEF *fitsbdf(mfd,mff,mfn,outflg)
/*+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
.PURPOSE   compute scaling information MIDAS image
.COMMENT   if mfd<0 the pointer is returned directly!    
.RETURN    return pointer to scaling structure, on error NULL pointer
---------------------------------------------------------------------*/
int        mfd;                /* IN: MIDAS file number              */
int        mff;                /* IN: MIDAS file data format         */
char       *mfn;               /* IN: name of MIDAS file             */
int        *outflg;	       /* OUT: 99, if BSCALE, BZERO are used
				       0, if LHCUTS was set or computed */
{
  int            i, n, ns, nb, nv;
  int            mmm, ioff, null, naxis;
  int            nax[MXDIM], unit[4];
  void           fpeh();
  float          *f, *fb, cuts[2];
  double         lcut, hcut;


  *outflg = 0;

  if (mfd<0) return &sdef;


  SCDRDI(mfd,"NAXIS",1,1,&nv,&naxis,unit,&null);
  sdef.dsize = (naxis) ? 1 : 0;
  SCDRDI(mfd,"NPIX",1,MXDIM,&nv,nax,unit,&null);
  for (i=0; i<naxis; i++) sdef.dsize *= nax[i];

  cuts[0] = cuts[1] = 0.0;
  SCDRDR(mfd,"LHCUTS",3,2,&nv,cuts,unit,&null);
  lcut = cuts[0]; hcut = cuts[1];

  if (nv!=2 || hcut<=lcut) {        	/* compute max/min cuts of data */

    struct FCT_STRUCT *fctpntr;
    int   istat;
    float  bscale, bzero;

    fctpntr = FCT.ENTRIES + mfd;
    if (fctpntr->FILTYP > 0) 
       {				/* take FITS header as it is */
       istat = SCDRDR(mfd,"BSCALE",1,1,&nv,&bscale,unit,&null);
       if (istat == ERR_NORMAL) 
          istat = SCDRDR(mfd,"BZERO",1,1,&nv,&bzero,unit,&null);
       if (istat == ERR_NORMAL) 
          {
          sdef.bfac = bscale;
          sdef.boff = bzero;
          }
       else
          {
          sdef.bfac = 1.0;
          sdef.boff = 0.0;
          }

       if (nv == 2)			/* low, hi cut are equal */
          {
          sdef.dmin = sdef.dmax = cuts[0];
          *outflg = 0;
          }
       else
          *outflg = 99;
       return &sdef;
       } 

    fb = (float *)osmmget(BSIZE);	/* get internal buffer */
    if (!fb) {
      SCTPUT("Error: cannot allocate internal buffer");
      SCFCLO(mfd); return (SDEF *) 0;
    }
    osscatch(SIGFPE,fpeh); fpexc = 0;       /* chatch fp-exception  */
    ioff = 1; nb = BSIZE/4;
    lcut = hcut = 0.0;
    if (mff!=D_R4_FORMAT) {
      SCFCLO(mfd);
      SCFOPN(mfn,D_R4_FORMAT,0,F_IMA_TYPE,&mmm);
    }
    else mmm = mfd;
    while (ioff<=sdef.dsize) {       /* find max/min values in data */
      f = fb;
      n = (sdef.dsize-ioff<nb) ? sdef.dsize-ioff+1 : nb;
      SCFGET(mmm, ioff, n, &ns, (char *) f);
      if (ioff == 1) lcut = hcut = *f;
      ioff += ns;
      if (ns<n) ioff = sdef.dsize + 1;    /* to force exit of loop  */
      while (ns--) {
	fpexc = 0;
	if (!isNULLFLOAT(*f))
	  if (*f<lcut && !fpexc) lcut = *f;
	  else if (hcut<*f && !fpexc) hcut = *f;
	f++;
      }
    }
    if (mff!=D_R4_FORMAT) {
      SCFCLO(mmm);
      SCFOPN(mfn,mff,0,F_IMA_TYPE,&mfd);
    }
    osmmfree((char *) fb);
  }

  sdef.dmax = hcut;
  sdef.dmin = lcut;
  sdef.bfac = 0.5 * (hcut-lcut)/2147483645.0;
  if (sdef.bfac==0.0) sdef.bfac = 1.0;
  sdef.boff = 0.5 * (hcut + lcut);

  return &sdef;
}

TXDEF *fitstbl(mfd,ffmt,cut)
/*+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
.PURPOSE   get table information
.COMMENT   if mfd<0 the pointer is returned directly!    
.RETURN    return pointer to Table structure
---------------------------------------------------------------------*/
int        mfd;                /* IN: MIDAS file number              */
int        ffmt;               /* IN: file format  B/O/X             */
int        cut;                /* IN: cut flag                       */
{
  char          *pc, fc;
  unsigned char  bnull;
  short          snull;
  int            ns, i, nn, dtype;
  int            len, fr, fw, fm;
  int            bytes, items, null;
  FDEF           *fd;

  if (mfd<0) return &txdef;

  TCIGET(mfd,&(txdef.tfields),&(txdef.nrow),&ns,&ns,&ns);
  txdef.mxrow = 0; txdef.mxcol = 0; txdef.col = fdef;
  fd = fdef;
  for (i=1; i<=txdef.tfields; i++, fd++) {   /* get info. for columns */
    TCFGET(mfd,i,fd->tdisp,&len,&dtype);
    TCBGET(mfd,i,&dtype,&items,&bytes);
    pc = fd->tdisp;
    while (*pc && *pc!=' ') pc++;
    *pc = '\0';
    if (ffmt == 'B') fd->trepn =  (dtype!=D_C_FORMAT) ? 1 : bytes;
    else fd->trepn =  (dtype!=D_C_FORMAT) ? items : bytes;  
    strcpy(fd->tform,fd->tdisp);
    dcffmt(fd->tform,&fr,&fc,&fw,&fm);
    TCLGET(mfd,i,fd->ttype);
    TCUGET(mfd,i,fd->tunit);
    fd->tnnul = 0;
    switch (dtype) {
       case D_R4_FORMAT :
            fd->tdfmt = 'E';
            if (ffmt!='B') fd->twdth = 4*items;
	    else if (cut=='C') {
	      if (fc=='I') fd->tdfmt = 'I';
	      fd->twdth = fw;
	    }
	    else { strcpy(fd->tform,"E15.8"); fd->twdth = 15; }
	    break;
       case D_R8_FORMAT :
            fd->tdfmt = 'D';
	    if (ffmt!='B') fd->twdth = 8*items;
	    else if (cut=='C') fd->twdth = fw;
	    else {
	      strcpy(fd->tform,"D24.16"); fd->twdth = 24;
	    }
	    break;
       case D_C_FORMAT :
	    fd->tdfmt = 'A';
	    if (ffmt!='B') fd->twdth = bytes;
	    else if (cut=='C') fd->twdth = fw;
	    else {
	      sprintf(fd->tform,"A%d",bytes); fd->twdth = bytes;
	    }
	    break;
       case D_I1_FORMAT :
	    fd->tdfmt = 'S';
	    if (ffmt!='B') fd->twdth = 2*items;
	    else if (cut=='C') fd->twdth = fw;
	    else {
	      strcpy(fd->tform,"I4"); fd->twdth = 4;
	    }
	    TBL_toNULL((TBL_D_I2<<TBL_D_BITS)|1, (char *) &snull);
	    fd->tnnul = snull;
	    break;
       case D_I2_FORMAT :
	    fd->tdfmt = 'S';
	    if (ffmt!='B') fd->twdth = 2*items;
	    else if (cut=='C') fd->twdth = fw;
	    else {
	      strcpy(fd->tform,"I6"); fd->twdth = 6;
	    }
	    TBL_toNULL((TBL_D_I2<<TBL_D_BITS)|1, (char *) &snull);
	    fd->tnnul = snull;
	    break;
       case D_I4_FORMAT :
	    fd->tdfmt = 'I';
	    if (ffmt!='B') fd->twdth = 4*items;
	    else if (cut=='C') fd->twdth = fw;
	    else {
	      strcpy(fd->tform,"I11"); fd->twdth = 11;
	    }
	    TBL_toNULL((TBL_D_I4<<TBL_D_BITS)|1, (char *) &null);
	    fd->tnnul = null;
	    break;
     }
    txdef.mxrow += fd->twdth;
    if (txdef.mxcol<fd->twdth) txdef.mxcol = fd->twdth;
  }

  return &txdef;
}

void fpeh(s)                                 /* fp-exception handler   */
int  s;
{
  osscatch(SIGFPE,fpeh);
  fpexc++;
}