C @(#)darithm.for	13.1.1.1 (ESO-DMD) 06/02/98 18:17:33
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 Massachusetts Ave, Cambridge, 
C MA 02139, USA.
C
C Correspondence 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
      PROGRAM DARITHM
C
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C.IDENTIFICATION
C  program DARITHM		version 1.00	850822
C  K. Banse                  	ESO - Garching
C
C.KEYWORDS
C  double precision, arithmetic operations, bulk data frames
C
C.PURPOSE
C  evaluate an arithmetic expression involving frames, constants and functions
C  and store the result into a data frame or key OUTPUTR(1), if used as pocket
C  calculator
C
C.ALGORITHM
C  "clean" the expression by replacing all frame names by F, all constants by C
C  and all functions by P, convert it to polish (postfix) notation and evaluate
C  it piecewise
C  establish special condition handler for arithmetic traps...
C  copy all descriptors of first frame operand to result frame.
C
C.INPUT/OUTPUT
C  the following keys are used:
C
C  DEFAULT/C/1/2	1. elem. = Y, result goes to frame
C			         = N, result is a constant -> key OUTPUTR(1)
C			2. elem. = I (use world coords),
C 				 = P (use frame pixels, i.e frames are arrays)
C 				 = X (use frame pixels + work on planes)
C  OUT_A/C/1/80		result frame, if DEFAULT(1:1) = Y
C  Pi/C/1/100		arithmetic expression, may contain up to
C			12 operands;
C			i = 3 or 1, if DEFAULT = Y or N
C	
C.VERSIONS
C  see SCCS
C	
C-------------------------------------------------------------------------
C
      IMPLICIT NONE
C
      INTEGER      NULCNT,TSTFLG
      INTEGER      FRACNT
      INTEGER      STAT,IAV,N,OFF
      INTEGER      UNIT(1),NULLO
      INTEGER      MADRID(1)
C
      REAL         RNULL(2)
C 
      DOUBLE PRECISION      USRNUL,DEPSP,DEPSN
C
      CHARACTER    ZLINE*320,ZRESFR*100,ZDEFA*2
      CHARACTER    PARMS(8)*2,CBUF*100
C 
      INCLUDE  'MID_INCLUDE:ST_DEF.INC'
C 
      COMMON      /NULCOM/ NULCNT,TSTFLG,USRNUL,DEPSP,DEPSN
      COMMON      /IPOOL/  FRACNT
C
      COMMON      /VMR/ MADRID
C
      INCLUDE  'MID_INCLUDE:ST_DAT.INC'
C 
      DATA   ZRESFR  /' '/
      DATA   ZLINE   /' '/
      DATA   PARMS   /'P1','P2','P3','P4','P5','P6','P7','P8'/
C
C  set up MIDAS environment
      CALL STSPRO('DARITHM ')
C
C  get default indicator + null data
      CALL STKRDC('DEFAULT',1,1,3,IAV,ZDEFA,UNIT,NULLO,STAT)
      CALL STKRDR('NULL',2,2,IAV,RNULL,UNIT,NULLO,STAT)
      NULCNT = 0         
      USRNUL = DBLE(RNULL(1))
      IF (RNULL(2).GT.0.) THEN
         TSTFLG = 1
      ELSE
         TSTFLG = 0
      ENDIF
      DEPSP = 10.E-49
      DEPSN = -DEPSP
C
C  get expression
      OFF = 1
      IF (ZDEFA(1:1).EQ.'N') THEN
         DO 400, N=1,8
            CALL STKRDC(PARMS(N),1,1,100,IAV,CBUF,UNIT,NULLO,STAT)
            IF (CBUF(1:1).EQ.'?') GOTO 1000
C 
            IF (N.EQ.8) THEN                  !P8 may have embedded blanks
               ZLINE(OFF:) = CBUF(1:)//' '
               GOTO 1010
            ENDIF
C 
            IAV = INDEX(CBUF,' ')             !no blanks in here ...
            IF (IAV.LE.1) IAV = 100
            ZLINE(OFF:) = CBUF(1:)
            OFF = OFF + IAV - 1
400      CONTINUE
C 
C  here for result = expression
C 
      ELSE                                               !get result frame...
         CALL STKRDC('OUT_A',1,1,100,IAV,ZRESFR,UNIT,NULLO,STAT)
         DO 800, N=3,8
            CALL STKRDC(PARMS(N),1,1,100,IAV,CBUF,UNIT,NULLO,STAT)
            IF (CBUF(1:1).EQ.'?') GOTO 1000
C 
            IF (N.EQ.8) THEN                  !P8 may have embedded blanks
               ZLINE(OFF:) = CBUF(1:)//' '
               GOTO 1010
            ENDIF
C
            IAV = INDEX(CBUF,' ')             !no blanks in here ...
            IF (IAV.LE.1) IAV = 100
            ZLINE(OFF:) = CBUF(1:)
            OFF = OFF + IAV - 1
800      CONTINUE
      ENDIF
C
1000  ZLINE(OFF:) = ' '
C
1010  CALL STKRDC('HISTORY',1,10,4,IAV,CBUF,UNIT,NULLO,STAT)  !get qualifier
      CALL UPCAS(CBUF(1:4),CBUF(1:4))
      IF ((CBUF(1:2).NE.'IM') .AND. (CBUF(1:2).NE.'PI')) THEN
         CALL STETER(47,'Command/qualifier currently NOT supported')
      ELSE
         CALL ARTHMD(ZDEFA,ZRESFR,ZLINE)
      ENDIF
C
      RNULL(1) = FLOAT(NULCNT)
      CALL STKWRR('NULL',RNULL,1,1,UNIT,STAT)     !update null count in key NULL
      IF (NULCNT.GT.0) CALL DDSPNL(NULCNT)        !display if NULL values there
C 
      CALL STSEPI
      END