C @(#)smoot.for 17.1.1.1 (ES0-DMD) 01/25/02 17:17:57 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 SUBROUTINE SMOOT (RMAT,NY,NX) C C 11.X.88 C IMPLICIT NONE INTEGER NX INTEGER NY REAL RMAT(NY,NX) INTEGER MX1 INTEGER MY1 PARAMETER (MX1=55) PARAMETER (MY1=55) C INTEGER N35, N36 INTEGER I, J, L, K, INC, INR REAL SP, S REAL TMT(3,3) REAL RMAT2(MY1,MX1) C DATA TMT/1.0,1.0,1.0,1.0,2.0,1.0,1.0,1.0,1.0/ C N35 = 3 N36 = N35/2+1 SP = 10. DO 100 K = N36,NY-N36+1 INC = K-N36 DO 110 I = N36,NX-N36+1 S = 0. INR = I-N36 DO 120 J = 1,N35 DO 130 L = 1,N35 S = S+RMAT(INR+L,INC+J)*TMT(L,J) 130 CONTINUE 120 CONTINUE RMAT2(I,K) = S/SP 110 CONTINUE 100 CONTINUE C DO 200 K = N36,NY-N36+1 DO 210 I = N36,NX-N36+1 RMAT(I,K) = RMAT2(I,K) 210 CONTINUE 200 CONTINUE RETURN END REAL FUNCTION GALE2(NPX,NPY,H,SIG,BET,XA,YA,X0,Y0) C++++++++++ C CALCOLA L' INTEGRALE BIDIMENSIONALE DI UNA C CURVA DI MOFFAT IN UN PIXEL CON IL METODO C DI GAUSS-LEGENDRE C C H = ALTEZZA CURVA DI MOFFAT C SIG = R " " " C BET = BETA " " " C X0 - Y0 = COORD. CENTRO " " " C XA - YA = COORD. INIZIALI DI INTEGRAZIONE C NPX = NUMERO DI PUNTI INTERMEDI SULL' ASSE X C SU CUI CALCOLARE L' INTEGRALE NUMERICO (MAX=6) C NPY = NUMERO DI PUNTI INTERMEDI SULL' ASSE Y C SU CUI CALCOLARE L' INTEGRALE NUMERICO (MAX=6) C------------- IMPLICIT NONE INTEGER NPX INTEGER NPY REAL H REAL SIG REAL BET REAL XA REAL YA REAL X0 REAL Y0 C REAL X, Y REAL XB, YB REAL XT, YT REAL VT REAL PSF INTEGER I, J REAL HE(6,6),HX(6,6) C DATA HX/ 0.0 ,0.0 ,0.0 ,0.0 ,0.0 ,0.0, * -.577350,.577350 ,0.0 ,0.0 ,0.0 ,0.0, * -.774597, 0. ,.774597 ,0.0 ,0.0 ,0.0, * -.861136,-.339981,.339981 ,.861136,0.0 ,0.0, * -.906179,-.538469,0. ,.538469,.906179,0.0, * -.932469,-.661209,-.238619,.238619,.661209,.932469/ DATA HE/2. ,0.0 ,0.0 ,0.0 ,0.0 ,0.0, * 1. ,1. ,0.0 ,0.0 ,0.0 ,0.0, * .555555,.888889,.555555,0.0 ,0.0 ,0.0, * .347855,.652145,.652145,.347855,0.0 ,0.0, * .236927,.478629,.568889,.478629,.236927,0.0, * .171324,.360762,.467914,.467914,.360762,.171324/ C PSF(H,X,Y,X0,Y0,SIG,BET) = H*(1+((X-X0)**2+(Y-Y0)**2) * /SIG**2)**(-BET) XB = XA+1 YB = YA+1 VT = 0. C DO 10 I = 1,MIN0(6,NPX) XT = .5*(HX(I,NPX)+XB+XA) DO 20 J = 1,MIN0(6,NPY) YT = .5*(HX(J,NPY)+YB+YA) VT = VT+HE(I,NPX)*HE(J,NPY)*PSF(H,XT,YT,X0,Y0,SIG,BET) 20 CONTINUE 10 CONTINUE GALE2 = VT/4. RETURN END REAL FUNCTION GALEG(NPX,NPY,H,SIG,XA,YA,X0,Y0) C C CALCOLA L' INTEGRALE BIDIMENSIONALE DI UNA C CURVA DI MOFFAT IN UN PIXEL CON IL METODO C DI GAUSS-LEGENDRE C C H = ALTEZZA CURVA DI MOFFAT C SIG = R " " " C X0 - Y0 = COORD. CENTRO " " " C XA - YA = COORD. INIZIALI DI INTEGRAZIONE C NPX = NUMERO DI PUNTI INTERMEDI SULL' ASSE X C SU CUI CALCOLARE L' INTEGRALE NUMERICO (MAX=6) C NPY = NUMERO DI PUNTI INTERMEDI SULL' ASSE Y C SU CUI CALCOLARE L' INTEGRALE NUMERICO (MAX=6) C IMPLICIT NONE INTEGER NPX INTEGER NPY REAL H REAL SIG REAL XA REAL YA REAL X0 REAL Y0 C REAL X, Y REAL XB, YB REAL XT, YT REAL VT REAL PSF INTEGER I, J REAL HE(6,6),HX(6,6) C DATA HX/ 0.0 ,0.0 ,0.0 ,0.0 ,0.0 ,0.0, * -.577350,.577350 ,0.0 ,0.0 ,0.0 ,0.0, * -.774597, 0. ,.774597 ,0.0 ,0.0 ,0.0, * -.861136,-.339981,.339981 ,.861136,0.0 ,0.0, * -.906179,-.538469,0. ,.538469,.906179,0.0, * -.932469,-.661209,-.238619,.238619,.661209,.932469/ DATA HE/2. ,0.0 ,0.0 ,0.0 ,0.0 ,0.0, * 1. ,1. ,0.0 ,0.0 ,0.0 ,0.0, * .555555,.888889,.555555,0.0 ,0.0 ,0.0, * .347855,.652145,.652145,.347855,0.0 ,0.0, * .236927,.478629,.568889,.478629,.236927,0.0, * .171324,.360762,.467914,.467914,.360762,.171324/ C PSF(H,X,Y,X0,Y0,SIG) = H*EXP(-4*ALOG(2.)*((X-X0)**2+(Y-Y0)**2) * /SIG**2) XB = XA+1 YB = YA+1 VT = 0. C DO 10 I=1,MIN0(6,NPX) XT = .5*(HX(I,NPX)+XB+XA) DO 20 J=1,MIN0(6,NPY) YT = .5*(HX(J,NPY)+YB+YA) VT = VT+HE(I,NPX)*HE(J,NPY)*PSF(H,XT,YT,X0,Y0,SIG) 20 CONTINUE 10 CONTINUE C GALEG=VT/4. RETURN END