C @(#)starsa.for	19.1 (ES0-DMD) 02/25/03 13:25:40
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+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C
C-----------------------------------------------------------------------
      SUBROUTINE STARSA(IOPC, AS, JAPYS, IBUFS, I,
     &                  J, LPXL, LSBP, NCT, PMT,
     &                  PRC, APSF, FPSF, NPAS, SCALE,
     &                  CINT)
C
      IMPLICIT  NONE
      INCLUDE   'MID_REL_INCL:INVENT.INC/NOLIST'
C
      INTEGER   IOPC               ! IN: -1 subtract , 1 add
      REAL      AS(1)              ! MOD: Data array
      INTEGER   JAPYS(1)           ! IN: Pointers to data lines
      INTEGER   IBUFS(4)           ! IN: Limits of data region
      INTEGER   I                  ! IN: Object x coordinate
      INTEGER   J                  ! IN: Object y coordinate
      INTEGER   LPXL               ! IN: Pixels extend
      INTEGER   LSBP               ! IN: Subpixels extend
      INTEGER   NCT(NIPAR)         ! IN: Component integer parameters
      REAL      PMT(NRPAR)         ! IN: Component real parameters
      REAL      PRC(0:MAXSUB)      ! IN: Object profile
      REAL      APSF(0:MAXSUB)     ! IN: One dimensional p.s.f.
      REAL      FPSF(1)            ! IN: Two dimensional p.s.f.
      INTEGER   NPAS               ! IN: Iteration count
      REAL      SCALE              ! IN: Central point scaling
      REAL      CINT               ! OUT: Component contribution
C
      INTEGER   IADR , IADR0 , IARG , IDIS , IK , IK1 , IK2 , IN , IS
      INTEGER   JADR , JK , JK1 , JK2 , JN , JOFF , JOFS , JS
      INTEGER   KSAT , LN , MPXL , MPXL2 , MSBP
C
      REAL   ADIS , AI , AJ , ATMP , CORR , DIST , SCA
C
      LOGICAL   CENTER , UNDER
C     LOGICAL   CENTR
C
      CINT = 0.0
      IN = NCT(1)
      JN = NCT(2)
      LN = NCT(5)
      KSAT = NCT(6)
      IK1 = MAX( IN-LN , IBUFS(1)+I )
      IK2 = MIN( IN+LN , IBUFS(3)+I )
      JK1 = MAX( JN-LN , IBUFS(2)+J )
      JK2 = MIN( JN+LN , IBUFS(4)+J )
      AI = PMT(10)
      AJ = PMT(11)
      IF ( NPAS .GT. 2 ) THEN
          SCA = PMT(12)
      ELSE
          IF ( KSAT .EQ. -1 ) THEN
              SCA = SCALE * PMT(2)
          ELSE
              ATMP = APSF(KSAT+2)
              IF ( ATMP .GT. 0.0 ) THEN
                  SCA = PRC(KSAT+2) / ATMP
              ELSE
                  SCA = 0.0
              ENDIF
          ENDIF
      ENDIF
C
      MPXL = 2 * LPXL + 1
      MSBP = 2 * LSBP + 1
      MPXL2 = MPXL * MPXL
C
      IF ( LPXL .GT. 0 .OR. LSBP .GT. 0 ) THEN
          UNDER = .TRUE.
          IS = NINT( ( AI - FLOAT(IN) ) * FLOAT(MSBP) )
          JS = NINT( ( AJ - FLOAT(JN) ) * FLOAT(MSBP) )
          IF ( IS .LT. (-LSBP) ) THEN
              IS = -LSBP
          ELSE IF ( IS .GT. LSBP ) THEN
              IS = LSBP
          ENDIF
          IF ( JS .LT. (-LSBP) ) THEN
              JS = LSBP
          ELSE IF ( JS .GT. LSBP ) THEN
              JS = LSBP
          ENDIF
          IADR0 = ( (LSBP+JS) * MSBP + LSBP + IS ) * MPXL2 +
     &                   LPXL * (MPXL+1) + 1
      ELSE
          UNDER = .FALSE.
      ENDIF
C
      JOFF = J + IBUFS(2) - 1
      DO 40 JK = JK1 , JK2
          IF ( UNDER .AND. IABS(JK-JN) .LE. LPXL ) THEN
              CENTER = .TRUE.
              JADR = IADR0 + (JK-JN)*MPXL
          ELSE
              CENTER = .FALSE.
          ENDIF
          JOFS = JAPYS(JK-JOFF) - I
          DO 50 IK = IK1 , IK2
              IF ( CENTER .AND. IABS(IK-IN) .LE. LPXL ) THEN
                  IADR = JADR + IK-IN
                  CORR = FPSF(IADR) * SCA
              ELSE
                  DIST = SQRT( (AI-FLOAT(IK))**2.0 +
     &                              (AJ-FLOAT(JK))**2.0 )
                  IDIS = INT( DIST )
                  ADIS = DIST - IDIS
                  CORR = ( APSF(IDIS)*(1.0-ADIS) +
     &                             APSF(IDIS+1)*ADIS ) * SCA
              ENDIF
              IARG = JOFS + IK
              IF ( IOPC .EQ. -1 ) THEN
                  AS(IARG) = AS(IARG) - CORR
              ELSE IF ( IOPC .EQ. 1 ) THEN
                  AS(IARG) = AS(IARG) + CORR
              ENDIF
              IF ( JK .EQ. J .AND. IK .EQ. I ) THEN
                  CINT = CORR
              ENDIF
   50     CONTINUE
   40 CONTINUE
C
      RETURN
C
      END
