! @(#)tfocus.prg	17.1.1.1 (ES0-DMD) 01/25/02 17:16:55
!===========================================================================
! 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 Massachusetss Ave, Cambridge, 
! MA 02139, USA.
!
! Corresponding 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
!===========================================================================
!
! @(#)tfocus.prg	17.1.1.1 (ESO-ASD) 01/25/02 17:16:55
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++
!.IDENT      tfocus.prg
!.PURPOSE    Sequence focus in the displayed image
!.USE        @c tfocus
!.AUTHOR     Cristian Levin
!.VERSION    1.0       01/06/92  Creation
!.VERSION    2.0       12/08/94  plotting enhanced (CL)
!.VERSION    2.1       10/01/95  implementation of parameters file (CL)
!.VERSION    2.2       14/09/95  add CREATE/GRAPHICS + remove CLEAR/CHAN OVER
!				 K.Banse
!.KEYWORDS
! ------------------------------------------------------
!
! PARAMETERS:
!    P1: start-x (default = 0.)
!    P2: step-x  (default = 1.)
!
define/parameter p1 0. N "Enter start-x : "
define/parameter p2 1. N "Enter step-x : "
!
define/param p8 + ? "Enter input_spec for center/gauss:"
if p8(1:1) .eq. "+" then
   define/local centin/c/1/6 CURSOR
else
   define/local centin/c/1/60 {p8}
endif
!
! initialization
define/local focx/d/1/1 0.0
define/local focy/d/1/1 0.0
define/local sigmax/d/1/1 0.0
define/local sigmay/d/1/1 0.0
define/local focus/d/1/1 0.0
define/local aux/d/1/1 0.0
define/local fox/d/1/1 0.0
define/local sign/d/1/1 0.0
create/table tfocus 15 0 null
!
! setup the graphics window
create/graph
set/graphics default
set/graphics pmode=1
!
! getting the coordinates...
center/gauss {centin} tfocus
if outputi(1) .eq. 0 then
    write/out
    write/out Invalid number of coordinates.
    return
endif
!
set/midas output=no ! we don't want output on the screen...
!
!
! fit-x computation
compute/table tfocus :step = SEQ
regression/polynomial tfocus :xfwhm :step 2
save/regression tfocus test
comp/regression tfocus :fitx = test(:step)
focx = {p1} + {p2}*({outputd(2)}/(-2.*{outputd(3)})) - {p2}
!
! fit-y computation
regression/polynomial tfocus :yfwhm :step 2
save/regression tfocus test
comp/regression tfocus :fity = test(:step)
focy = {p1} + {p2}*({outputd(2)}/(-2.*{outputd(3)})) - {p2}
!
! compute the scaling for plots (big job ;{ )
!
compute/table tfocus :maxim = MAX(:xfwhm,:yfwhm)
sort/table tfocus :maxim(-)
define/local maxim/r/1/1 {tfocus,:maxim,1}
maxim = {maxim} + 0.2 	! let's give a margin value
compute/table tfocus :minim = MIN(:xfwhm,:yfwhm)
sort/table tfocus :minim
define/local minim/r/1/1 {tfocus,:minim,1}
minim = {minim} - 0.2 	! let's give a margin value
set/graphics yaxis={minim},{maxim}
!
! plot (x,y)-value versus (x,y)-sigma
sort/table tfocus :step
set/graphics stype=6 colour=1 ! cross(x), black
plot/table tfocus :step :xfwhm
set/graphics stype=5 colour=4 ! cross(+), blue
overplot/table tfocus :step :yfwhm
!
! plot of the fit
set/graphics stype=0 ltype=1 ! no symbol, solid line
set/graphics colour=1 ! black
overplot/table tfocus :step :fitx
set/graphics colour=4 ! blue
overplot/table tfocus :step :fity 
!
! focus computation
aux  = {focx} * {focy}
aux  = M$ABS({aux})
fox  = M$ABS({focx})
sign = {focx}/{fox}
focus = {sign} * M$SQRT({aux})
!
! residual computation
compute/table tfocus :residx2 = (:xfwhm - :fitx) ** 2
compute/table tfocus :residy2 = (:yfwhm - :fity) ** 2
statis/table tfocus :residx2
sigmax = outputr(3)
statis/table tfocus :residy2
sigmay = outputr(3)
!
set/midas output=yes ! enable the screen output...
write/out
write/out
if mid$sess(16:20) .eq. "94NOV" then
   set/format f12.2
else
   set/format ,f12.2			!for double keywords
endif
write/out "X: BLACK    y: BLUE"
write/out
write/out "FOCUS:  X =  {focx} " "Y =  {focy}  " "MEAN =  {focus}"
write/out 
if mid$sess(16:20) .eq. "94NOV" then
   set/format f12.3
else
   set/format ,f12.3			!for double keywords
endif
write/out "sigma-x = {sigmax}   sigma-y = {sigmay}"
set/graphics default
! 
outputr(1) = focx
outputr(2) = focy
outputr(3) = focus
outputr(4) = sigmax
outputr(5) = sigmay