C @(#)optop.for 17.1.1.1 (ES0-DMD) 01/25/02 17:55:13 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 OPTOP C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C.IDENTIFICATION: program OPTOP C.PURPOSE: to provide X and Y coordinates for plates to be used at the C Cassegrain focus of the 3.6mt telescope. C.AUTHOR: A. Gemmo Padova Department of Astronomy C.VERSION: 861101 G. Lund Creation C 900625 S. Cristiani Last update on HP machine C 910301 A. Gemmo New portable version running in MIDAS env. C 910524 A. Gemmo Modified for OPTOPUS context. C 910717 A. Gemmo Modified to handle "-00" in coordinates C of the center. DEC read as a string. C------------------------------------------------------------------------------ IMPLICIT NONE C C *** C INTEGER MADRID(1) INTEGER NRINP INTEGER NCINP INTEGER NROUT INTEGER NCOUT PARAMETER (NCINP=7) PARAMETER (NCOUT=9) PARAMETER (NROUT=300) INTEGER ISTAT,IAC,IAV INTEGER LFIEL1,LFIEL2,LFIEL3,LFIEL4 INTEGER LFIEL5,LFIEL6,LFIEL7 INTEGER TYP1,TYP2,TYP3,TYP4 INTEGER TYP5,TYP6,TYP7 INTEGER KUN,KNUL INTEGER NACINP,NARINP,NSINP INTEGER OUTTYP,OUTCOL(10),COLNUM(10),CN(10) INTEGER TIDINP INTEGER TIDOUT INTEGER NBIG INTEGER I1,I7,I8,I9,I10,NCOL,II INTEGER I12,I13,I14,I15,I16,I17,I18,J4,J5,J7 INTEGER K1,K2,NCINP1,IND1,IND2,IND3 C C *** C INTEGER ICOUN0,ICOUN1,ICOUN2,ICOUN3 INTEGER FLAG(NROUT),FLAG1(NROUT),FLAG2(NROUT) C C *** C REAL CDEC(3) C C *** C DOUBLE PRECISION FRAD,SGSEP,BGSEP,OOSEP,OGSEP DOUBLE PRECISION SCALFA,EXDGUI,EXDOBJ DOUBLE PRECISION XFAREA,YFAREA DOUBLE PRECISION CRA(3),DCDEC(3),CRA1,CDEC1 DOUBLE PRECISION DRA(NROUT),DDEC(NROUT) DOUBLE PRECISION X(NROUT),Y(NROUT),Z(NROUT) DOUBLE PRECISION PI DOUBLE PRECISION RARAD,DECRAD,CRARAD,CDECRA DOUBLE PRECISION DISTAN DOUBLE PRECISION FAR1,FAR2 DOUBLE PRECISION DBGUID,DSGUID DOUBLE PRECISION XX(360),YY(360) DOUBLE PRECISION RAXRAD(360),DECYRAD(360) DOUBLE PRECISION RAX(360),DECY(360) DOUBLE PRECISION XF(4),YF(4),XFRAD(4),YFRAD(4) DOUBLE PRECISION FX(4),FY(4) DOUBLE PRECISION EQUI0,EQUI1,PCRARA,PCDECR DOUBLE PRECISION PCRA,PCDEC C C *** CHARACTER*80 STRING(25) CHARACTER*60 INPFIL CHARACTER*60 OUTFIL CHARACTER*16 ILAB(NCINP),RLAB(NCINP),OUTLAB CHARACTER*16 OUTUNI CHARACTER*16 OUTFOR CHARACTER*8 FORM CHARACTER*16 FORMC16,FORMC1,FORMR8 CHARACTER*16 IDENT(NROUT) CHARACTER*1 TYPE(NROUT) CHARACTER*1 CHECK(NROUT) CHARACTER*80 CDECCH CHARACTER*80 DECCHA(3) CHARACTER*1 SIGN CHARACTER*1 ACFLAG,PFLAG C C *** LOGICAL NULL C C *** INCLUDE 'MID_INCLUDE:TABLES.INC' COMMON /VMR/ MADRID INCLUDE 'MID_INCLUDE:TABLED.INC' C C *** DATA ILAB/'IDENT ','TYPE ','CHECK ', 2 'OLDRA ','OLDDEC ','RA ','DEC '/ C C *** C DATA FORMC16/'A16'/ DATA FORMC1 /'A1'/ DATA FORMR8 /'F7.0'/ DATA PI/3.14159265358979D0/ DATA ICOUN0,ICOUN1,ICOUN2,ICOUN3/0,0,0,0/ DATA DBGUID,DSGUID/0.0D0,0.0D0/ DATA FX/-20010.,-20010.,20010.,20010./ DATA FY/137000.,118100.,118100.,137000./ DATA SIGN/'+'/ C C *** start definition of constraints of starplates C DATA FRAD /1.37D5/ !radius of the field (microns) DATA SGSEP /3.45D3/ !min. sep. between small guidestars (microns) DATA BGSEP /1.25D4/ !min. sep. between big guidestars (microns) DATA NBIG /2/ !number of big guidestars DATA OOSEP /3.45D3/ !min. object-object separation (microns) DATA OGSEP /9.00D3/ !min. object-guidestar separation (microns) DATA SCALFA /1.56D-7/ !scal. fact. for depth of the holes (microns) DATA EXDGUI /-6.00D2/ !extra-depth for guidestars (microns) DATA EXDOBJ /-1.49D3/ !extra-depth for objects (microns) DATA XFAREA /2.001D4/ !half X-dimension of forbidden area (microns) DATA YFAREA /1.890D4/ !Y-dimension of forbidden area (microns) C C *** start the code C CALL STSPRO('OPTOP') C DO II=1,NROUT FLAG (II) = 0 FLAG1(II) = 0 FLAG2(II) = 0 ENDDO C C *** get the input table and the output table C CALL STKRDC('INPUTFI',1,1,60,IAC,INPFIL,KUN,KNUL,ISTAT) CALL STKRDC('OUTPUTF',1,1,60,IAC,OUTFIL,KUN,KNUL,ISTAT) C C *** open input table C CALL TBTOPN(INPFIL,F_I_MODE,TIDINP,ISTAT) C C *** get information about the input table C CALL TBIGET(TIDINP,NCOL,NRINP,NSINP,NACINP,NARINP,ISTAT) IF(NRINP.EQ.0)THEN !no rows in table STRING(1) = '*** FATAL: There are no data in the input TABLE' CALL STETER(9,STRING(1)) ENDIF C CALL TBLSER(TIDINP,'IDENT',COLNUM(1),ISTAT) CALL TBLSER(TIDINP,'TYPE',COLNUM(2),ISTAT) CALL TBLSER(TIDINP,'CHECK',COLNUM(3),ISTAT) CALL TBLSER(TIDINP,'OLDRA',COLNUM(4),ISTAT) CALL TBLSER(TIDINP,'OLDDEC',COLNUM(5),ISTAT) CALL TBLSER(TIDINP,'RA',COLNUM(6),ISTAT) CALL TBLSER(TIDINP,'DEC',COLNUM(7),ISTAT) C DO I1=1,NCINP IF(COLNUM(I1).GE.1)THEN CALL TBLGET(TIDINP,COLNUM(I1),RLAB(I1),ISTAT) ELSE GOTO 1902 ENDIF IF(RLAB(I1).NE.ILAB(I1))THEN STRING(2) = '*** FATAL: Wrong column label' CALL STETER(9,STRING(2)) ENDIF 1902 ENDDO C IF (COLNUM(1).GE.1)THEN CALL TBFGET(TIDINP,COLNUM(1),FORM,LFIEL1,TYP1,ISTAT) IF (TYP1.NE.D_C_FORMAT)THEN STRING(3) = '*** FATAL: Wrong `IDENT` column format' CALL STETER(9,STRING(3)) ENDIF ENDIF IF (COLNUM(2).GE.1)THEN CALL TBFGET(TIDINP,COLNUM(2),FORM,LFIEL2,TYP2,ISTAT) IF (TYP2.NE.D_C_FORMAT)THEN STRING(4) = '*** FATAL: Wrong `TYPE` column format' CALL STETER(9,STRING(4)) ENDIF ENDIF IF (COLNUM(3).GE.1)THEN CALL TBFGET(TIDINP,COLNUM(3),FORM,LFIEL3,TYP3,ISTAT) IF (TYP3.NE.D_C_FORMAT)THEN STRING(5) = '*** FATAL: Wrong `CHECK` column format' CALL STETER(9,STRING(5)) ENDIF ENDIF IF (COLNUM(4).GE.1)THEN CALL TBFGET(TIDINP,COLNUM(4),FORM,LFIEL4,TYP4,ISTAT) IF (TYP4.NE.D_R8_FORMAT)THEN STRING(6) = '*** FATAL: Wrong `OLDRA` column format' CALL STETER(9,STRING(6)) ENDIF ENDIF IF (COLNUM(5).GE.1)THEN CALL TBFGET(TIDINP,COLNUM(5),FORM,LFIEL5,TYP5,ISTAT) IF (TYP5.NE.D_R8_FORMAT)THEN STRING(7) = '*** FATAL: Wrong `OLDDEC` column format' CALL STETER(9,STRING(7)) ENDIF ENDIF IF (COLNUM(6).GE.1)THEN CALL TBFGET(TIDINP,COLNUM(6),FORM,LFIEL6,TYP6,ISTAT) IF (TYP6.NE.D_R8_FORMAT)THEN STRING(8) = '*** FATAL: `RA` column format' CALL STETER(9,STRING(8)) ENDIF ENDIF IF (COLNUM(7).GE.1)THEN CALL TBFGET(TIDINP,COLNUM(7),FORM,LFIEL7,TYP7,ISTAT) IF (TYP7.NE.D_R8_FORMAT)THEN STRING(9) = '*** FATAL: Wrong `DEC` column format' CALL STETER(9,STRING(9)) ENDIF ENDIF C C *** create the output table C CALL TBTINI(OUTFIL,0,F_O_MODE,NCOUT,NROUT,TIDOUT,ISTAT) C OUTTYP = TYP1 !initialize ident column OUTFOR = FORMC16 OUTUNI = ' ' OUTLAB = 'IDENT ' CALL TBCINI(TIDOUT,OUTTYP,16,OUTFOR,OUTUNI, 2 OUTLAB,OUTCOL(1),ISTAT) C OUTTYP = D_C_FORMAT !initialize type column OUTFOR = FORMC1 OUTUNI = ' ' OUTLAB = 'TYPE ' CALL TBCINI(TIDOUT,OUTTYP,1,OUTFOR,OUTUNI, 2 OUTLAB,OUTCOL(2),ISTAT) C OUTTYP = D_C_FORMAT !initialize check column OUTFOR = FORMC1 OUTUNI = ' ' OUTLAB = 'CHECK ' CALL TBCINI(TIDOUT,OUTTYP,1,OUTFOR,OUTUNI, 2 OUTLAB,OUTCOL(3),ISTAT) C OUTTYP = D_R8_FORMAT !initialize x column OUTFOR = FORMR8 OUTUNI = 'MICRONS' OUTLAB = 'X ' CALL TBCINI(TIDOUT,OUTTYP,1,OUTFOR,OUTUNI, 2 OUTLAB,OUTCOL(4),ISTAT) C OUTTYP = D_R8_FORMAT !initialize y column OUTFOR = FORMR8 OUTUNI = 'MICRONS' OUTLAB = 'Y ' CALL TBCINI(TIDOUT,OUTTYP,1,OUTFOR,OUTUNI, 2 OUTLAB,OUTCOL(5),ISTAT) C OUTTYP = D_R8_FORMAT !initialize z column OUTFOR = FORMR8 OUTUNI = 'MICRONS' OUTLAB = 'Z ' CALL TBCINI(TIDOUT,OUTTYP,1,OUTFOR,OUTUNI, 2 OUTLAB,OUTCOL(6),ISTAT) C C *** read center coordinates (if present) C CALL STKRDD('PLATECEN',1,1,IAV,CRA(1),KUN,KNUL,ISTAT) CALL STKRDD('PLATECEN',2,1,IAV,CRA(2),KUN,KNUL,ISTAT) CALL STKRDD('PLATECEN',3,1,IAV,CRA(3),KUN,KNUL,ISTAT) CALL STKRDC('ACFLAG',1,1,1,IAC,ACFLAG,KUN,KNUL,ISTAT) IF(ACFLAG.EQ.'Y'.OR.ACFLAG.EQ.'y')THEN CALL STKRDD('PLATECEN',4,1,IAV,DCDEC(1),KUN,KNUL,ISTAT) CALL STKRDD('PLATECEN',5,1,IAV,DCDEC(2),KUN,KNUL,ISTAT) CALL STKRDD('PLATECEN',6,1,IAV,DCDEC(3),KUN,KNUL,ISTAT) GOTO 8888 ELSE IF(ACFLAG.EQ.'N'.OR.ACFLAG.EQ.'n')THEN CALL STKRDC('PLATECHA',1,1,80,IAV,CDECCH,KUN,KNUL,ISTAT) ENDIF C IND1 = INDEX(CDECCH,',') IND2 = INDEX(CDECCH(IND1+1:80),',') + IND1 DECCHA(1) = CDECCH(1:IND1-1) IND3 = INDEX(CDECCH,'-') DECCHA(2) = CDECCH(IND1+1:IND2-1) DECCHA(3) = CDECCH(IND2+1:80) C C *** decode character DECCHA to real array C CALL USRINP(CDEC(1),80,'R',DECCHA(1)) CALL USRINP(CDEC(2),80,'R',DECCHA(2)) CALL USRINP(CDEC(3),80,'R',DECCHA(3)) C C *** write decoded dec in keyword C DCDEC(1) = DBLE(CDEC(1)) DCDEC(2) = DBLE(CDEC(2)) DCDEC(3) = DBLE(CDEC(3)) C IF(DCDEC(1).LT.-90..OR.DCDEC(1).GT.90.)THEN CALL STETER(ISTAT,'*** FATAL: Wrong DEC-degrees !!') ELSE IF(DCDEC(2).GT.60.)THEN CALL STETER(ISTAT,'*** FATAL: Wrong DEC-arcmin !!') ELSE IF(DCDEC(3).GT.60.)THEN CALL STETER(ISTAT,'*** FATAL: Wrong DEC-arcsec !!') ENDIF C CALL STKWRD('PLATECEN',DCDEC(1),4,1,KUN,ISTAT) CALL STKWRD('PLATECEN',DCDEC(2),5,1,KUN,ISTAT) CALL STKWRD('PLATECEN',DCDEC(3),6,1,KUN,ISTAT) C C *** convert RA and DEC to decimal format C 8888 CRA1 = CRA(1)+CRA(2)/60.+CRA(3)/3600. IF(DCDEC(1).LT.0..OR. 1 DCDEC(2).LT.0..OR. 2 DCDEC(3).LT.0..OR. 3 IND3.NE.0) THEN SIGN(1:1) = '-' DCDEC(1) = DABS(DCDEC(1)) DCDEC(2) = DABS(DCDEC(2)) DCDEC(3) = DABS(DCDEC(3)) ENDIF CDEC1 = DCDEC(1)+DCDEC(2)/60.+DCDEC(3)/3600. IF(SIGN(1:1).EQ.'-')CDEC1 = -CDEC1 C C *** do the work C C C *** start computation of X and Y positions C DO I7=1,NRINP !read in table elements CALL TBERDC(TIDINP,I7,COLNUM(1),IDENT(I7),KNUL,ISTAT) CALL TBERDC(TIDINP,I7,COLNUM(2),TYPE(I7),KNUL,ISTAT) C C *** convert TYPE to uppercase C IF(TYPE(I7).EQ.'o')TYPE(I7) = 'O' IF(TYPE(I7).EQ.'s')TYPE(I7) = 'S' IF(TYPE(I7).EQ.'b')TYPE(I7) = 'B' C CALL TBERDD(TIDINP,I7,COLNUM(6),DRA(I7),KNUL,ISTAT) CALL TBERDD(TIDINP,I7,COLNUM(7),DDEC(I7),KNUL,ISTAT) ENDDO C C *** convert center RA and DEC to radians C CRARAD = CRA1*1.5D1*(PI/1.8D2) CDECRA = CDEC1*(PI/1.8D2) C C *** precess center coordinates C CALL STKRDC('PFLAG',1,1,1,IAC,PFLAG,KUN,KNUL,ISTAT) CALL STKRDD('EQUINOX',1,1,IAC,EQUI0,KUN,KNUL,ISTAT) CALL STKRDD('EQUINOX',2,1,IAC,EQUI1,KUN,KNUL,ISTAT) C IF(PFLAG.EQ.'N'.OR.PFLAG.EQ.'n')THEN CALL STTPUT(' ',ISTAT) CALL STTPUT('No precession applied to center coordinates', 1ISTAT) CALL STTPUT(' ',ISTAT) PCRARA = CRARAD PCDECR = CDECRA ELSE IF(PFLAG.EQ.'Y'.OR.PFLAG.EQ.'y')THEN CALL PRE(CRARAD,CDECRA,PCRARA,PCDECR,EQUI0,EQUI1) ENDIF C C *** convert precessed RA and DEC to hours and degrees C PCRA = PCRARA*(1.8D2/PI)/1.5D1 PCDEC = PCDECR*(1.8D2/PI) CALL STKWRD('PLATEC1',PCRA,1,1,KUN,ISTAT) CALL STKWRD('PLATEC1',PCDEC,2,1,KUN,ISTAT) C C *** subroutine ADXY taken from Yanling Fang's astrometric package C (modified by A. Gemmo) DO J7=1,NRINP RARAD = DRA(J7)*1.5D1*(PI/1.8D2) DECRAD = DDEC(J7)*(PI/1.8D2) CALL ADXY(RARAD,DECRAD,PCRARA,PCDECR,X(J7),Y(J7)) ENDDO C C *** fill in table elements C DO J4=1,NRINP CALL TBEWRC(TIDOUT,J4,OUTCOL(1),IDENT(J4),ISTAT) CALL TBEWRC(TIDOUT,J4,OUTCOL(2),TYPE(J4),ISTAT) CALL TBEWRD(TIDOUT,J4,OUTCOL(4),X(J4),ISTAT) CALL TBEWRD(TIDOUT,J4,OUTCOL(5),Y(J4),ISTAT) ENDDO C C *** check maximum distance from the center C DO I8=1,NRINP DISTAN = DSQRT(X(I8)**2+Y(I8)**2) IF(DISTAN.GE.FRAD)THEN ICOUN1 = ICOUN1+1 WRITE(STRING(7),222)IDENT(I8) WRITE(STRING(8),333)TYPE(I8),X(I8),Y(I8) 222 FORMAT(1X,'Object: ',A16) 333 FORMAT(1X,'Type: ',A1,',',' X = ',F10.0,' Y = ',F10.0, 2 ' is too far from the center') CALL STTPUT(STRING(7),ISTAT) CALL STTPUT(STRING(8),ISTAT) CHECK(I8) = 'N' ENDIF ENDDO IF(ICOUN1.NE.0)THEN WRITE(STRING(9),444)ICOUN1,NRINP 444 FORMAT(1X,'You are going to loose',I3,' object(s) out of ',I3) CALL STTPUT(STRING(9),ISTAT) STRING(10)=' You can delete it/them or change the center positi 2on' CALL STTPUT(STRING(10),ISTAT) ENDIF C C *** check if there is any object in the forbidden area C FAR1 = XFAREA FAR2 = FRAD-YFAREA DO I9=1,NRINP IF(DABS(X(I9)).LE.FAR1.AND.Y(I9).GE.FAR2)THEN ICOUN2 = ICOUN2+1 WRITE(STRING(11),555)IDENT(I9) WRITE(STRING(12),666)TYPE(I9),X(I9),Y(I9) 555 FORMAT(1X,'Object: ',A16) 666 FORMAT(1X,'Type: ',A1,',',' X = ',F10.0,' Y = ',F10.0, 2 ' is in the forbidden area') CALL STTPUT(STRING(11),ISTAT) CALL STTPUT(STRING(12),ISTAT) CHECK(I9) = 'N' ENDIF ENDDO IF(ICOUN2.NE.0)THEN WRITE(STRING(13),777)ICOUN2 777 FORMAT(1X,I3,' object(s) is/are in the forbidden area') CALL STTPUT(STRING(13),ISTAT) WRITE(STRING(13),778) 778 FORMAT(1X,' You MUST delete it/them or change the center posit 2ion') CALL STTPUT(STRING(13),ISTAT) ENDIF C C *** check the number of big guidestars C DO I10=1,NRINP IF(TYPE(I10).EQ.'B')THEN ICOUN3 = ICOUN3+1 ENDIF ENDDO IF(ICOUN3.LT.NBIG)THEN STRING(14)= ' *** FATAL: Your `BIG` guidestars are not enough' CALL STETER(9,STRING(14)) ENDIF C C *** check the reciprocal superpositions of big guidestars C DO I12=1,NRINP IF(TYPE(I12).EQ.'B'.AND.FLAG(I12).EQ.0)THEN DO I13=1,NRINP IF(TYPE(I13).EQ.'B'.AND.I12.NE.I13.AND. 2 FLAG(I13).EQ.0)THEN DISTAN = DSQRT((X(I12)-X(I13))**2 2 +(Y(I12)-Y(I13))**2) DBGUID = DMAX1(DBGUID,DISTAN) IF(DISTAN.LT.BGSEP)THEN WRITE(STRING(16),888)IDENT(I12),IDENT(I13) 888 FORMAT(1X,' Big guidestars: ',A16,' and ',A16,' are too close') CALL STTPUT(STRING(16),ISTAT) CHECK(I12) = 'N' CHECK(I13) = 'N' FLAG(I12) = FLAG(I12)+1 FLAG(I13) = FLAG(I13)+1 ENDIF ELSE GOTO 1111 ENDIF 1111 ENDDO ELSE GOTO 1112 ENDIF 1112 ENDDO C C *** check the reciprocal superpositions of small guidestars and objects C DO I14=1,NRINP IF(TYPE(I14).NE.'B'.AND.FLAG1(I14).EQ.0)THEN DO I15=1,NRINP IF(TYPE(I15).NE.'B'.AND.I14.NE.I15.AND. 2 FLAG1(I15).EQ.0)THEN DISTAN = DSQRT((X(I14)-X(I15))**2 2 +(Y(I14)-Y(I15))**2) DSGUID = DMAX1(DSGUID,DISTAN) IF (DISTAN.LT.SGSEP)THEN IF(TYPE(I14).EQ.'S'.AND.TYPE(I15).EQ.'S')THEN WRITE(STRING(18),2222)IDENT(I14),IDENT(I15) 2222 FORMAT(1X,' Small guidestars: ',A16,' and ',A16, 2 ' are too close') CALL STTPUT(STRING(18),ISTAT) CHECK(I14) = 'N' CHECK(I15) = 'N' FLAG1(I14) = FLAG1(I14)+1 FLAG1(I15) = FLAG1(I15)+1 ELSE IF(TYPE(I14).EQ.'S'.AND.TYPE(I15).EQ.'O')THEN WRITE(STRING(18),2223)IDENT(I14),IDENT(I15) 2223 FORMAT(1X,' Small guidestar: ',A16,' and object:', A16) CALL STTPUT(STRING(18),ISTAT) WRITE(STRING(18),3333) 3333 FORMAT(1X,' are too close') CALL STTPUT(STRING(18),ISTAT) CHECK(I14) = 'N' CHECK(I15) = 'N' FLAG1(I14) = FLAG1(I14)+1 FLAG1(I15) = FLAG1(I15)+1 ELSE IF(TYPE(I14).EQ.'O'.AND.TYPE(I15).EQ.'S')THEN WRITE(STRING(18),2224)IDENT(I14),IDENT(I15) 2224 FORMAT(1X,' Object: ',A16,' and small guidestar:', A16) CALL STTPUT(STRING(18),ISTAT) WRITE(STRING(18),3334) 3334 FORMAT(1X,' are too close') CALL STTPUT(STRING(18),ISTAT) CHECK(I14) = 'N' CHECK(I15) = 'N' FLAG1(I14) = FLAG1(I14)+1 FLAG1(I15) = FLAG1(I15)+1 ELSE IF(TYPE(I14).EQ.'O'.AND.TYPE(I15).EQ.'O')THEN WRITE(STRING(18),2225)IDENT(I14),IDENT(I15) 2225 FORMAT(1X,' Objects: ',A16,' and ',A16,' are too close') CALL STTPUT(STRING(18),ISTAT) CHECK(I14) = 'N' CHECK(I15) = 'N' FLAG1(I14) = FLAG1(I14)+1 FLAG1(I15) = FLAG1(I15)+1 ENDIF ENDIF ELSE GOTO 1113 ENDIF 1113 ENDDO ELSE GOTO 1114 ENDIF 1114 ENDDO C C *** check the reciprocal superpositions of small and big guidestars C DO I16=1,NRINP IF(TYPE(I16).EQ.'B'.AND.FLAG2(I16).EQ.0)THEN DO I17=1,NRINP IF(TYPE(I17).NE.'B'.AND.I16.NE.I17.AND. 2 FLAG2(I17).EQ.0)THEN DISTAN = DSQRT((X(I16)-X(I17))**2 2 +(Y(I16)-Y(I17))**2) DSGUID = DMAX1(DSGUID,DISTAN) IF(DISTAN.LT.OGSEP)THEN IF (TYPE(I16).EQ.'B'.AND.TYPE(I17).EQ.'S')THEN WRITE(STRING(20),4444)IDENT(I16),IDENT(I17) 4444 FORMAT(1X,' Big guidestar: ',A16, ' and small guidestar: ',A16) CALL STTPUT(STRING(20),ISTAT) WRITE(STRING(21),5555) 5555 FORMAT(1X,' are too close') CALL STTPUT(STRING(21),ISTAT) CHECK(I16) = 'N' CHECK(I17) = 'N' ELSE IF(TYPE(I16).EQ.'B'.AND.TYPE(I17).EQ.'O')THEN WRITE(STRING(20),4441)IDENT(I16),IDENT(I17) 4441 FORMAT(1X,' Big guidestar: ',A16,' and object: ',A16) WRITE(STRING(21),5551) 5551 FORMAT(1X,' are too close') CALL STTPUT(STRING(21),ISTAT) CHECK(I16) = 'N' CHECK(I17) = 'N' FLAG2(I16) = FLAG2(I16)+1 FLAG2(I17) = FLAG2(I17)+1 ENDIF ENDIF ELSE GOTO 4442 ENDIF 4442 ENDDO ELSE GOTO 4443 ENDIF 4443 ENDDO C C *** compute depth of the hole for big guidestars, small guidestars C and objects C DO I18=1,NRINP IF(TYPE(I18).EQ.'B')THEN Z(I18) = -(X(I18)*X(I18)+Y(I18)*Y(I18))*SCALFA+EXDGUI ELSE IF(TYPE(I18).NE.'B')THEN Z(I18) = -(X(I18)*X(I18)+Y(I18)*Y(I18))*SCALFA+EXDOBJ ENDIF X(I18)=-X(I18) ENDDO C C *** fill in the output table C DO J5=1,NRINP CALL TBEWRC(TIDOUT,J5,OUTCOL(3),CHECK(J5),ISTAT) CALL TBEWRD(TIDOUT,J5,OUTCOL(6),Z(J5),ISTAT) ENDDO C C *** calculate circle to show plate limits C CALL TBTOPN('circle',F_IO_MODE,TIDINP,ISTAT) !open table C CALL TBIGET(TIDINP,NCINP1,NRINP,NSINP,NACINP,NARINP,ISTAT) !table info C DO K2=1,NRINP !fill in RAX and DECY columns in circle.tbl C C *** subroutine XYAD taken from Yanling Fang's astrometric package C (God help us!) CALL TBLSER(TIDINP,'XX',CN(1),ISTAT) CALL TBLSER(TIDINP,'YY',CN(2),ISTAT) CALL TBERDD(TIDINP,K2,CN(1),XX(K2),NULL,ISTAT) CALL TBERDD(TIDINP,K2,CN(2),YY(K2),NULL,ISTAT) CALL XYAD(XX(K2),YY(K2),PCRARA,PCDECR, 2 RAXRAD(K2),DECYRAD(K2)) RAX(K2) = RAXRAD(K2)*360./(2.*15.*PI) DECY(K2) = DECYRAD(K2)*360./(2.*PI) CALL TBLSER(TIDINP,'RAX',CN(3),ISTAT) CALL TBLSER(TIDINP,'DECY',CN(4),ISTAT) CALL TBEWRD(TIDINP,K2,CN(3),RAX(K2),ISTAT) CALL TBEWRD(TIDINP,K2,CN(4),DECY(K2),ISTAT) ENDDO C CALL STKWRD('CIRDIM',RAX(1),1,1,KUN,ISTAT) CALL STKWRD('CIRDIM',DECY(90),2,1,KUN,ISTAT) C C *** conversion of forbidden area coordinates C DO K1=1,4 CALL XYAD(FX(K1),FY(K1),PCRARA,PCDECR,XFRAD(K1),YFRAD(K1)) XF(K1) = XFRAD(K1)*360./(2.*15.*PI) CALL STKWRD('XF',XF(K1),K1,1,KUN,ISTAT) YF(K1) = YFRAD(K1)*360./(2.*PI) CALL STKWRD('YF',YF(K1),K1,1,KUN,ISTAT) ENDDO C C *** initialize row selection flag C CALL TBSINI(TIDOUT,ISTAT) C C *** close the table C CALL TBTCLO(TIDOUT,ISTAT) C C *** over and out C CALL STSEPI END C C *** Subroutines used in the code: C SUBROUTINE PRE(A0,D0,A1,D1,EQX0,EQX1) C C ... calculate general precession for two given epochs C C ... References: STUMPFF P.,ASTRON. ASTROPHYS. SUPPL. 41, P. 1 (1980) C ... EXPLAN. SUPPL. ASTRON. EPHEREMERIS, P. 28 (1961) C IMPLICIT NONE !REAL*8 (A-H,O-Z) C DOUBLE PRECISION A0,D0,D1,EQX0,EQX1 DOUBLE PRECISION DT0,DT,DTS,DTC,DC1900,DC1M2 DOUBLE PRECISION DC1,DC2,DC3,DC4,DC5,DC6,DC7,DC8,DC9 DOUBLE PRECISION DZED,DZETA,DCSAR,DTHETA DOUBLE PRECISION PI,SINR,COSR,A1,R C DATA DCSAR/4.848136812D-6/,DC1900/1900.0D0/,DC1M2/0.01D0/, * DC1/2304.25D0/,DC2/1.396D0/,DC3/0.302D0/,DC4/0.018D0/, * DC5/0.791D0/,DC6/2004.683D0/,DC7/-0.853D0/,DC8/-0.426D0/, * DC9/-0.042D0/ C PI = 3.1415926535897928 C C ... first determine the precession angles zeta, z and theta C DT0=(EQX0-DC1900)*DC1M2 DT=(EQX1-EQX0)*DC1M2 DTS=DT*DT DTC=DTS*DT DZETA=((DC1+DC2*DT0)*DT+DC3*DTS+DC4*DTC)*DCSAR DZED=DZETA + DC5*DTS*DCSAR DTHETA=((DC6+DC7*DT0)*DT+DC8*DTS+DC9*DTC)*DCSAR C C ... with these and A0,D0 compute the precessed coordinates A1,D1 C D1 = ASIN(DCOS(DTHETA)*DSIN(D0)+DSIN(DTHETA)*DCOS(D0)* 1 DCOS(A0+DZETA)) SINR = (DCOS(D0)*DSIN(A0+DZETA))/DCOS(D1) COSR = (DCOS(DTHETA)*DCOS(D0)*DCOS(A0+DZETA)- 1 DSIN(DTHETA)*DSIN(D0))/DCOS(D1) R = DASIN(SINR) IF (COSR.LT.0) R=PI-R ! determine quadrant of R A1 = R+DZED ! this is the precessed R.A. IF (A1.GT.2.0D0*PI) A1=A1-2.0D0*PI IF (A1.LT.0.0D0) A1=A1+2.0D0*PI C RETURN END C C *************** XYAD ******* C C The subroutine XYAD transforms X and Y coordinates in RAs and DECs C SUBROUTINE XYAD(XX,YY,RA0RAD,D0RAD,RARAD,DECRAD) C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C XX (INPUT) Value of X coordinate to be transformed (microns) C YY (INPUT) Value of Y coordinate to be transformed (microns) C RA0RAD (INPUT) Center RA (radians) C D0RAD (INPUT) Center DEC (radians) C RARAD (OUTPUT) RA (radians) C DECRAD (OUTPUT) DEC (radians) C------------------------------------------------------------------------------ IMPLICIT NONE C C *** C DOUBLE PRECISION PI,FOC DOUBLE PRECISION XX,YY,AKSI,ANU DOUBLE PRECISION COSD0,SIND0 DOUBLE PRECISION CX,CY,BX(9),BY(9) DOUBLE PRECISION TERM1,TERM2,XX1,YY1 DOUBLE PRECISION RH,SINPH,COSPH,SINRH,COSRH DOUBLE PRECISION RARAD,DECRAD DOUBLE PRECISION RA0RAD,D0RAD C C *** C DATA PI/3.14159265358979D0/ DATA BX/0.0D0,0.0D0,0.0D0,0.0D0,0.0D0, 2 0.0D0,0.0D0,0.0D0,0.0D0/ DATA BY/0.0D0,0.0D0,0.0D0,0.0D0,0.0D0, 2 0.0D0,0.0D0,0.0D0,0.0D0/ DATA CX,CY/0.0D0,0.0D0/ C C *** do the work C FOC = 1.0D0/3.6D3 BX(1) = 9.9973D-1/0.14006 BY(2) = BX(1) C COSD0 = DCOS(D0RAD) SIND0 = DSIN(D0RAD) C TERM1 = BX(1)*BY(2)-BX(2)*BY(1) TERM2 = BX(2)*BY(1)-BX(1)*BY(2) C XX1 = XX*FOC/1.0D3 YY1 = YY*FOC/1.0D3 C ANU = XX1*BY(1)-(YY1*BY(2)*TERM2/TERM1) 2 -(CY*BX(2)*BY(1)/TERM1)+(CY*BX(1)*BY(2)/TERM1) AKSI = (XX1*TERM1+ANU*BX(2)-CY*BX(2))/BY(2)+CX C RH = DSQRT(AKSI**2 + ANU**2) C SINPH = AKSI/RH COSPH = ANU/RH C RH = RH*PI/1.8D2 C SINRH= DSIN(RH) COSRH = DCOS(RH) C DECRAD = COSRH*SIND0+SINRH*COSD0*COSPH DECRAD = DATAN(DECRAD/DSQRT(1.0D0-DECRAD**2)) C RARAD = SINRH*SINPH/DCOS(DECRAD) RARAD = RA0RAD + DATAN(RARAD/DSQRT(1.0D0-RARAD**2)) C RETURN END C C C *************** ADXY ******* C C The subroutine ADXY transforms RAs and DECs in X and Y coordinates C SUBROUTINE ADXY(ALFAS,DELTS,ALFA0,DELT0,XX,YY) C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C ALFAS (INPUT) Value of RA to be transformed (radians) C DELTS (INPUT) Value of DEC to be transformed (radians) C ALFA0 (INPUT) RA of plate center (radians) C DELT0 (INPUT) DEC of plate center (radians) C XX (OUTPUT) X coordinate (microns) C YY (OUTPUT) Y coordinate (microns) C-------------------------------------------------------------------------- IMPLICIT NONE C C *** C INTEGER IDX(9),IDY(9) INTEGER NX,NY,IXY,NITER INTEGER I16,I17,N,L C C *** C DOUBLE PRECISION XX,YY DOUBLE PRECISION DELTS,ALFAS DOUBLE PRECISION RH,COSRH,SINRH,COSPH,SINPH DOUBLE PRECISION V,AKSI1,ANU1,AKSI,ANU DOUBLE PRECISION DKSI,DNU,DDIST,DCRIT DOUBLE PRECISION DKDX,DKDY,DNDX,DNDY DOUBLE PRECISION XBIS(10),AR(50,4),CX,CY,BX(9),BY(9) DOUBLE PRECISION ALFA0,DELT0,COSD0,SIND0 DOUBLE PRECISION XEST(2,2),YEST(2,2) DOUBLE PRECISION PI,FOC C C *** peculiar Cassegrain focus parameters C DATA NX/3/ DATA NY/3/ DATA CX/0.0D0/ DATA CY/0.0D0/ DATA IDX/1,1,1,0,0,0,0,0,0/ DATA IDY/1,1,1,0,0,0,0,0,0/ DATA BX/0.0D0,0.0D0,0.0D0,0.0D0,0.0D0, 2 0.0D0,0.0D0,0.0D0,0.0D0/ DATA BY/0.0D0,0.0D0,0.0D0,0.0D0,0.0D0, 2 0.0D0,0.0D0,0.0D0,0.0D0/ C C *** C DATA PI/3.14159265358979D0/ C C *** do the work C NITER = 0 FOC = 1.0D0/3.6D3 BX(1) = 9.9973D-1/0.14006 BY(2) = BX(1) C COSD0 = DCOS(DELT0) SIND0 = DSIN(DELT0) C IF (ALFA0.LT.PI/2D0.AND.ALFAS.GT.1.5D0*PI) ALFAS = ALFAS-2D0*PI IF (ALFA0.GT.1.5D0*PI.AND.ALFAS.LT.PI/2D0) ALFAS = ALFAS+2D0*PI C COSRH = COSD0*DCOS(DELTS)*DCOS(ALFAS-ALFA0) + SIND0*DSIN(DELTS) SINRH = DSQRT(1.0D0-COSRH**2) SINPH = DCOS(DELTS)*DSIN(ALFAS-ALFA0)/SINRH COSPH = (DSIN(DELTS)-SIND0*COSRH)/(COSD0*SINRH) RH = DATAN(SINRH/COSRH) C C *** estimate X and Y position C AKSI = RH*SINPH/(PI/1.8D2) ANU = RH*COSPH/(PI/1.8D2) C XEST(1,1) = ((AKSI-CX)*BY(2)-(ANU-CY)*BX(2))/ 2 (BX(1)*BY(2)-BX(2)*BY(1)) YEST(1,1) = ((AKSI-CX)*BY(1)-(ANU-CY)*BX(1))/ 2 (BX(2)*BY(1)-BX(1)*BY(2)) 30 XEST(2,1) = XEST(1,1)*1.0D3/FOC YEST(2,1) = YEST(1,1)*1.0D3/FOC C C60 NITER = NITER + 1 NITER = NITER + 1 C C70 AR(1,1) = XEST(1,1) AR(1,1) = XEST(1,1) AR(1,2) = YEST(1,1) AR(1,3) = 0.0D0 AR(1,4) = 0.0D0 C I16 = 1 N = NX+1 IXY = 1 CALL RAR(XBIS,N,IXY,IDX,IDY,AR,I16,0) L = NX V = 0.0D0 DO I17 = 1,L V = V + BX(I17)*XBIS(I17) ENDDO AKSI1 = V + CX C N = NY + 1 IXY = 2 CALL RAR(XBIS,N,IXY,IDX,IDY,AR,I16,0) L = NY V = 0.0D0 DO I17 = 1,L V = V + BY(I17)*XBIS(I17) ENDDO ANU1 = V + CY DKSI = AKSI1 - AKSI DNU = ANU1 - ANU DDIST = DSQRT(DKSI**2 + DNU**2) DCRIT = 1.0D-6 IF (DDIST.LE.DCRIT ) GOTO 1000 C C100 N = NX + 1 N = NX + 1 IXY = 1 CALL RAR(XBIS,N,IXY,IDX,IDY,AR,1,1) L = NX V = 0.0D0 DO I17 = 1,L V = V + BX(I17)*XBIS(I17) ENDDO DKDX = V C N = NX + 1 IXY = 1 CALL RAR(XBIS,N,IXY,IDX,IDY,AR,1,2) L = NX V = 0.0D0 DO I17 = 1,L V = V + BX(I17)*XBIS(I17) ENDDO DKDY = V C N = NY + 1 IXY = 2 CALL RAR(XBIS,N,IXY,IDX,IDY,AR,1,1) L = NY V = 0.0D0 DO I17 = 1,L V = V + BY(I17)*XBIS(I17) ENDDO DNDX = V C N = NY + 1 IXY = 2 CALL RAR(XBIS,N,IXY,IDX,IDY,AR,1,2) L = NY V = 0.0D0 DO I17 = 1,L V = V + BY(I17)*XBIS(I17) ENDDO DNDY = V C YEST(1,2) = (DKSI*DNDX-DNU*DKDX)/(DNDX*DKDY-DKDX*DNDY) XEST(1,2) = (DKSI*DNDY-DNU*DKDY)/(DKDX*DNDY-DNDX*DKDY) XEST(1,1) = XEST(1,1) - XEST(1,2) YEST(1,1) = YEST(1,1) - YEST(1,2) GOTO 30 1000 XX = -XEST(2,1) YY = YEST(2,1) RETURN END C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C C The subroutine RAR comes directly from the POS package (cfr Astrom. prgs.) C C---------------------------------------------------------------------------- SUBROUTINE RAR(XBIS,N,IXY,IDX,IDY,AR,I,K) C C *** C IMPLICIT NONE C C *** C INTEGER IDX(1),IDY(1) INTEGER K,I,I1,IXY,J,N C C *** C DOUBLE PRECISION X1(10) DOUBLE PRECISION XBIS(10),AR(50,4) C C *** C CHARACTER STRING*80 C C *** do the work C IF (K.EQ.0) GOTO 1 IF (K.EQ.1) GOTO 2 IF (K.EQ.2) GOTO 3 C 1 X1(1) = AR(I,1) X1(2) = AR(I,2) X1(3) = AR(I,1)*AR(I,2) X1(4) = AR(I,1)**2 X1(5) = AR(I,2)**2 X1(6) = AR(I,1)**3 X1(7) = AR(I,2)**3 X1(8) = AR(I,1)*AR(I,2)**2 X1(9) = AR(I,2)*AR(I,1)**2 GOTO 5 C 2 X1(1) = 1.0D0 X1(2) = 0.0D0 X1(3) = AR(I,2) X1(4) = 2.0D0*AR(I,1) X1(5) = 0.0D0 X1(6) = 3.0D0*AR(I,1)**2 X1(7) = 0.0D0 X1(8) = AR(I,2)**2 X1(9) = 2.0D0*AR(I,2)*AR(I,1) GOTO 5 C 3 X1(1) = 0.0D0 X1(2) = 1.0D0 X1(3) = AR(I,1) X1(4) = 0.0D0 X1(5) = 2.0D0*AR(I,2) X1(6) = 0.0D0 X1(7) = 3.0D0*AR(I,2)**2 X1(8) = 2.0D0*AR(I,1)*AR(I,2) X1(9) = AR(I,1)**2 C 5 IF(IXY.EQ.1) X1(10) = AR(I,3) IF(IXY.EQ.2) X1(10) = AR(I,4) IF(IXY.EQ.2) GOTO 20 I1 = 0 DO J=1,9 IF(IDX(J).EQ.0) GOTO 10 I1 = I1+1 XBIS(I1) = X1(J) 10 ENDDO IF (N-I1.NE.1) GOTO 50 XBIS(I1+1) = X1(10) GOTO 40 20 I1= 0 DO J=1,9 IF (IDY(J).EQ.0) GOTO 30 I1 = I1+1 XBIS(I1) = X1(J) 30 ENDDO IF(N-I1.NE.1) GOTO 50 XBIS(I1+1) = X1(10) 40 RETURN 50 STRING = '*** FATAL: Coeff. trouble' CALL STETER(9,STRING) RETURN END C---------------------------------------------------------------------------