C @(#)tsabnd.for	17.1.1.1 (ES0-DMD) 01/25/02 17:20:36
C===========================================================================
C Copyright (C) 1995 European Southern Observatory (ESO)
C
C This program is free software; you can redistribute it and/or 
C modify it under the terms of the GNU General Public License as 
C published by the Free Software Foundation; either version 2 of 
C the License, or (at your option) any later version.
C
C This program is distributed in the hope that it will be useful,
C but WITHOUT ANY WARRANTY; without even the implied warranty of
C MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
C GNU General Public License for more details.
C
C You should have received a copy of the GNU General Public 
C License along with this program; if not, write to the Free 
C Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, 
C MA 02139, USA.
C
C Corresponding concerning ESO-MIDAS should be addressed as follows:
C	Internet e-mail: midas@eso.org
C	Postal address: European Southern Observatory
C			Data Management Division 
C			Karl-Schwarzschild-Strasse 2
C			D 85748 Garching bei Muenchen 
C			GERMANY
C===========================================================================
C
C  @(#)tsabnd.for	5.1 (ESO-IPG) 4/5/93 15:58:48
C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C.COPYRIGHT (c)  1992     European Southern Observatory and Warsaw Observatory
C.IDENT     tsabnd.for
C.AUTHOR    Alex Schwarzenberg-Czerny, ESO and Warsaw Observatory
C.KEYWORD   MIDAS, time series analysis, BAND/TIME
C.LANGUAGE  FORTRAN 77
C.PURPOSE   Evaluate suitable frequency band for time analysis
C.VERSION   0.0               June 1992
C.RETURNS   None
C.ENVIRON   TSA context
C-----------------------------------------------------------------------------
C
C
      INCLUDE 'MID_REL_INCL:TSA_DEF.INC'
      INCLUDE 'MID_INCLUDE:ST_DEF.INC'
C
      CHARACTER*60 INAME                !  NAME OF INPUT TABLE
      INTEGER      MAXOBS               !  MAXIMUM NUMBER OF SCANNED OBSERV.
C
      REAL*8       START,STEP           !  OUTPUT FREQUENCY GRID
      INTEGER      NSTEPS
C
      INTEGER      IACTS,KUN,KNUL
      INTEGER      TID,IWORK,ITIME
      INTEGER      NCOL,ICOL,NROW,IROW,ISOR
      INTEGER      LFIELD,TTYP
      INTEGER      PTIME,PWORK
      CHARACTER*10 FORM
      CHARACTER*80 TEXT
C
      INCLUDE 'MID_REL_INCL:TSA_DAT.INC'
      INCLUDE 'MID_INCLUDE:ST_DAT.INC'
C
C
C   Get parameters
C
      CALL STSPRO ('tsabnd')
      CALL STKRDC ('IN_A',   1,1,60,IACTS,INAME ,KUN,KNUL,ISTAT)
      CALL STKRDI ('MAXOBS',   1, 1,IACTS,MAXOBS,KUN,KNUL,ISTAT)
C
C   Map input data
C
      CALL TBTOPN (INAME,F_I_MODE,TID,ISTAT)
      CALL TBIGET (TID,NCOL,NROW,ISOR,ICOL,IROW,ISTAT)
      CALL TBLSER (TID,'TIME',ITIME,ISTAT)
      IF (ITIME.LT.0) THEN
        CALL STETER(3,'Column :TIME not found')
      ENDIF
      CALL TBFGET (TID,ITIME,FORM,LFIELD,TTYP,ISTAT)
      CALL TBDGET (TID,ISTORE,ISTAT)
      IF (ISTORE.NE.F_TRANS) THEN
        TEXT='Input table '//INAME//' stored not transposed'
        CALL STETER(1,TEXT)
      ENDIF
      IF (TTYP.NE.D_R8_FORMAT) THEN
        CALL STETER(2,'Column :TIME must be in DOUBLE PRECISION')
      ENDIF
      CALL TBCMAP (TID,ITIME,PTIME,ISTAT)
C
C   Create work table
C
      CALL STFCRE('ZZMIDWORK',D_R8_FORMAT,F_X_MODE,F_IMA_TYPE,
     $    NROW,IWORK,ISTAT)
      CALL STFMAP(IWORK,F_X_MODE,1,NROW,IACTS,PWORK,ISTAT)
C
C   Evaluate frequency band
C
      CALL TIMBAND(MADRID(PTIME),MADRID(PWORK),NROW,MAXOBS,
     $    START,STEP,NSTEPS)
C
C   Return keyword values
C
      CALL STKWRD('STARTTSA', START, 1,1,KUN,ISTAT)
      CALL STKWRD('STEPTSA',  STEP,  1,1,KUN,ISTAT)
      CALL STKWRI('NSTEPS',NSTEPS,1,1,KUN,ISTAT)
      CALL STTPUT(
     $   'Keywords STARTTSA,STEPTSA and NSTEPS are set now.'
     $   ,ISTAT)
C
C   Wind-up
C
      CALL STSEPI
C
      END
C
C
C
      SUBROUTINE TIMBAND(TIM,DEL,NOBS,MAXOBS,START,STEP,NSTEPS)
      INCLUDE 'MID_REL_INCL:TSA_DEF.INC'
      INTEGER      NOBS,NSTEPS,MAXOBS
      DOUBLE PRECISION       TIM(NOBS),DEL(NOBS)
      REAL*8       START,STEP
C
      INTEGER      NOBS1,IOBS,MAXSTEPS
      REAL*8       FINISH,STEPS,SIZE
      CHARACTER*80 TEXT
      INCLUDE 'MID_REL_INCL:TSA_CONST.INC'
      DATA         MAXSTEPS/30000/
C
      IF (MAXOBS.EQ.0) MAXOBS=NOBS
      NOBS1=MIN(NOBS-1,MAXOBS)
      IF (NOBS1.LT.5) THEN
        CALL STETER(10,'Too few observations or MAXOBS too small')
      ENDIF
C
      DO 1 IOBS=1,NOBS1
        DEL(IOBS)=TIM(IOBS+1)-TIM(IOBS)
        IF (DEL(IOBS).LT.ZERO) THEN
          CALL STETER(14,':TIME must be sorted in ascending order')
        ENDIF
 1    CONTINUE
      CALL SORT(NOBS1,DEL)
C
      STEP=TIM(NOBS)-TIM(1)
      IF (STEP.LE.ZERO) THEN
        CALL STETER(11,'Input table has wrong :TIME numbers')
      ENDIF
      STEP=0.3/STEP
C
      SIZE=LOG(ONE*NOBS*ONE)
      IOBS=NOBS1*(ONE/(6.+0.3*SIZE)+0.05)+1
      FINISH=DEL(IOBS)
      IF (FINISH.LE.ZERO) THEN
        CALL STETER(12,
     $ 'Too finely spaced observations: bin them coarsly')
      ENDIF
      FINISH=HALF/FINISH*(DEL(NOBS1/2)/FINISH)**(0.6)
C
      STEPS=FINISH/STEP
C
      CALL STTPUT(' RESULTS OF FREQUENCY BAND EVALUATION:',ISTAT)
      WRITE(TEXT,'(2(A,1PE10.1))') 'Max. Frequency: ',FINISH,
     $   '  Resolution: ',STEP
      CALL STTPUT(TEXT,ISTAT)
      WRITE(TEXT,'(A,1PE10.1)') 'No. of points:  ',STEPS
      CALL STTPUT(TEXT,ISTAT)
      IF (STEPS.GT.MAXSTEPS) THEN
        STEPS=MAXSTEPS
        TEXT='*** DANGER *** Data span too long interval'//
     $     ' for good sampling of periodogrammes.'
        CALL STTPUT(TEXT,ISTAT)
        TEXT='Analysing data split into shorter'//
     $        'intervals and taking'//
     $       ' average of periodogrammes '//
     $       'will help by reducing resolution.'
        CALL STTPUT(TEXT,ISTAT)
      ENDIF
C
      NSTEPS=STEPS
      STEP=FINISH/STEPS
      START=ZERO
      END
C
C
C