C @(#)rfotmodel.for 17.1.1.1 (ES0-DMD) 01/25/02 17:18:17 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 PROGRAM MODEL C+++ C.PURPOSE: Interpolation program for undersampling problem in ROMAFOT C.AUTHOR: R. Buonnanno C.REVISED: R.H. Warmels C.VERSION: C--------------- IMPLICIT NONE C INTEGER MTG(55,55) DOUBLE PRECISION DX,DY DOUBLE PRECISION IN,ER DOUBLE PRECISION DX0,DY0 DOUBLE PRECISION VOLU,VOLG,TNOIS INTEGER IAV, ISTAT INTEGER IDX, IDY INTEGER I, J, II INTEGER INCR INTEGER KUN, KNUL INTEGER KVX, KVY INTEGER KVM1, KVM2, KVM3, KVM INTEGER MC, MR INTEGER NCHAR INTEGER NP3, NPUG, NSA INTEGER NPX, NPY INTEGER NTOT, NPXS, NPYS C REAL AMD, AML REAL BAF REAL BET, DITE, EFF, FCON REAL COSTMA REAL ERROLD, ERR REAL FX, FY, FLUSSO, FSP REAL H REAL PMA, PPO, PG, PSF REAL QS REAL RINPUT(15) REAL RMA, RON, RNOIS REAL SEEING, SIG, SKYL REAL TESP, TET REAL XB, YB REAL X0, Y0 REAL XA, YA REAL VM, VX, VY, VT, VF REAL VFCG, VFCGS REAL GALE2 C CHARACTER FRE*1,STRING*80,FILE*60 C DATA PG/3.14159/ DATA SKYL/3.18E-06/ C 999 FORMAT('*** INFO: Ratio FWHM over Pixel size =', E11.6) 998 FORMAT(2X,I3,2X,I3,2X,A1) 997 FORMAT(2X,F5.3,3(2X,F6.2),2X,F9.4) 996 FORMAT(2(2X,F3.1),2X,F6.2) 995 FORMAT(25I3) C C *** See in MIDAS a alive CALL STSPRO('MODEL') C C *** the default setting: C DEFAULT: C SKYL=3.18E-06 photon surface brightness of the sky in C photons cm**-2 s**-1 A**-1 arcsec**-2 C corresponding to V=21.3 at 5500 A C SEEING1 = 1 seeing (in arcsec) C BET1 = 3 beta C DITE1 = 2.2 diameter of the telescope (metres) C BAF1 = 800. bandwidth C EFF1 = .3 total efficiency C FCON1 = 10 photons-to-ADU conversion factor C TESP1 = 1000. exposure time C TET1 = .35 pixel size (in arcsec) C RON = 35 read out noise of the detector C AMD1 = 23. limiting magnitude (faint end) C AML1 = 18. limiting magnitude (bright end) C PMA1 = .5 sampling in magnitude C PPO = .1 positional accuracy?? C C *** get the input parameters CALL STKRDR('INPUTR',1,13,IAV,RINPUT,KUN,KNUL,ISTAT) SEEING = RINPUT(1) ! seeing BET = RINPUT(2) DITE = RINPUT(3)*100. BAF = RINPUT(4) EFF = RINPUT(5) FCON = RINPUT(6) TESP = RINPUT(7) TET = RINPUT(8) RON = RINPUT(9) AMD = RINPUT(10) AML = RINPUT(11) PMA = RINPUT(12) QS = RINPUT(13) C FSP = SEEING/TET SIG = (SEEING/TET)/(2*SQRT(2**(1./BET)-1)) PPO = .1 COSTMA = -2.5*ALOG10(FCON/(DITE**2*BAF*TESP*EFF*PG/4.))+7.5 MR = 7+INT(FSP/1.4+.7)*2 IF (MR.GT.15) THEN MR = 15 ENDIF MC = MR IDX = 0 IDY = 0 X0 = MC/2 Y0 = MR/2 FRE = 'N' XA = INT(X0)-MC/2+IDX YA = INT(Y0)-MR/2+IDY XB = XA+MC YB = YA+MR C C *** WRITE(STRING,999) FSP CALL STTPUT(STRING,ISTAT) CALL STKRDC('INPUTC',1,1,60,IAV,FILE,KUN,KNUL,ISTAT) NCHAR = INDEX(FILE,' ')-1 FILE = FILE(1:NCHAR)//'.dat' OPEN(UNIT=1,STATUS='NEW',ACCESS='SEQUENTIAL',FORM='FORMATTED', 2 FILE=FILE) C C *** write the output table WRITE(1,998) MR,MC,FRE WRITE(1,997) PPO,AMD,AML,PMA,COSTMA C INCR = -INT(PPO*100+0.5) DO 1000 KVX = 50,0,INCR VX = FLOAT(KVX)/100. X0 = MC/2+VX INCR = -INT(PPO*100+0.5) DO 1100 KVY = 50,KVX,INCR VY = FLOAT(KVY)/100. Y0 = MR/2+VY KVM1 = INT(100*AML+0.5) KVM2 = INT(100*AMD+0.5) KVM3 = INT(100*PMA+0.5) DO 1200 KVM = KVM1,KVM2,KVM3 VM = FLOAT(KVM)/100. RMA = VM FLUSSO = 10.**(-(RMA-7.5)/2.5) H = (BET-1)*FLUSSO*DITE**2*BAF*TESP*EFF/ 2 (4*SIG**2*FCON) FX = X0-INT(X0) FY = Y0-INT(Y0) WRITE(1,996) FX,FY,RMA C C *** calcuate matrix NP3 = 0 C NPUG = NPU VOLU = 0. VOLG = 0. TNOIS = 0. DO 1300 I=1,MR YA=INT(Y0)-MR/2+I-1+IDY DO 1400 J=1,MC XA=INT(X0)-MC/2+J-1+IDX C C *** calcualte int vero IF (BET.EQ.1.5) THEN DX = XA+1-X0 DY = YA+1-Y0 DX0 = DX-1. DY0 = DY-1. ER = SIG VT = H*(IN(DX,DY,ER)-IN(DX,DY0,ER)- 2 IN(DX0,DY,ER)+IN(DX0,DY0,ER))*ER**2 ELSE VT = GALE2(4,4,H,SIG,BET,XA,YA,X0,Y0) END IF FLUSSO = VT*FCON/(EFF*DITE**2*BAF*TESP*PG/4.) RNOIS = SQRT((FLUSSO+SKYL*TET**2.)*EFF*TESP* 2 BAF*DITE**2.*PG/4+RON**2) C TNOIS = TNOIS+RNOIS**2 VOLU = VOLU+VT VF = PSF(H,XA+.5,YA+.5,X0,Y0,SIG,BET) C IF (FRE.EQ.'N') THEN VFCG = VF ERR = RNOIS*QS/FCON NPUG = 1 1401 IF (ABS(VT-VFCG).LE.ERR.OR.NPUG.GT.5) GO TO 1402 NPUG = NPUG+1 VFCG = GALE2(NPUG,NPUG,H,SIG,BET,XA,YA,X0,Y0) GO TO 1401 1402 CONTINUE VOLG = VOLG+VFCG MTG(I,J) = NPUG**2 NP3 = NP3+NPUG**2 ELSE ERR = RNOIS*QS/FCON NSA = 1000 ERROLD = 10.E15 DO 1500 NPX=1,4 DO 1600 NPY=1,4 VFCG=GALE2(NPX,NPY,H,SIG,BET,XA,YA,X0,Y0) IF (ABS(VT-VFCG).LT.ERR) THEN NTOT=NPX*NPY IF (NTOT.LE.NSA) THEN IF (ERR.LT.ERROLD. OR. 2 NTOT.NE.NSA) THEN NSA = NTOT NPXS = NPX NPYS = NPY VFCGS = VFCG ERROLD = ERR END IF END IF END IF 1600 CONTINUE 1500 CONTINUE C VFCG = VFCGS VOLG = VOLG+VFCG MTG(I,J) = 10*NPXS+NPYS NP3 = NP3+NSA END IF 1400 CONTINUE 1300 CONTINUE C DO 1201 II=1,MR WRITE(1,995) (MTG(II,J),J=1,MC) 1201 CONTINUE 1200 CONTINUE 1100 CONTINUE 1000 CONTINUE C CLOSE(UNIT=1) CALL STSEPI END REAL FUNCTION PSF(H,X,Y,X0,Y0,SIG,BET) C REAL H REAL X REAL Y REAL X0 REAL Y0 REAL SIG REAL BET C PSF = H * (1+((X-X0)**2 + (Y-Y0)**2)/SIG**2)**(-BET) C RETURN END DOUBLE PRECISION FUNCTION IN(DX,DY,ER) C DOUBLE PRECISION DX DOUBLE PRECISION DY DOUBLE PRECISION ER C IN = DATAN(DX*DY/(ER*DSQRT(DX**2+DY**2+ER**2))) C RETURN END