C
        PROGRAM VES
C
C     Para compilar:
C     f77 -o ves ves.f -L$PGPLOT_DIR -lpgplot -lX11
C----------------------------------------------------------------------
C
       IMPLICIT NONE
C
       INTEGER I,J,N,NN
       INTEGER NWIN          ! Numero de ventanas activas en el plot.
       PARAMETER (NWIN=19)
       INTEGER NOBS          ! Cuantos observatorios
       PARAMETER (NOBS=10)
       INTEGER NNOBS         ! Observatorio elegido
       INTEGER PGOPEN
       INTEGER INDATE(3),INTIME(3)
       INTEGER IDAY,IMONTH,IYEAR
       INTEGER IHOUR,IMIN,ISEC
       INTEGER NSTAR,NSTAROLD  ! Numero de objetos trazados.
       INTEGER TESTN  
       INTEGER DEV(2)
C
       REAL WIN(4,NWIN)
       REAL ANGLE            ! Valores para PGPTXT
       REAL XX,YY            ! Posicion del cursor
       REAL X,Y
       REAL X1,X2,Y1,Y2
       REAL*8 TS0              ! Tiempo sidereo en Greenwich a 0h TU
       REAL*8 TS1              ! TS0 al dia siguiente.
       REAL*8 FJ0            ! Fecha Juliana a O TU 
       REAL*8 FJ1            ! FJ0+1. Fecha Juliana de "maniana" 
       REAL*8 FJ             ! Fecha Juliana para dia IDAY y hora RHOUR
       REAL*8 ARSUN,DECSUN     ! AR Y DEC del Sol en la fech de observ.(FJ)
       REAL*8 SEMIA,SEMIB      ! semiejes del elipsoide WGS-84
       REAL*8 REARTH           ! radio terrestre para el observatorio
       REAL*8 ANGHOA           ! angulo depresion horizonte debido a la altura
       REAL*8 PI
       PARAMETER(PI=3.141592654D0)
       REAL*8 DISZEN           ! distancia zenital del sol
       REAL*8 LAT,LONG,ALTOBS  ! Latitud, longitud y altura del obser.(DD.dddd)
       REAL*8 HORTO,HOCASO     ! angulos horarios de ortos y ocasos
       REAL*8 AORTO,AOCASO     ! azimuts de ortos y ocasos
       REAL*8 TUORTO1,TUORTO2     ! Tiempo universal de orto y ocaso
       REAL*8 TUOCASO1,TUOCASO2   ! utilizados para la iteracion
       REAL*8 TSL,TSL0             ! Tiempo sidereo local, y a 0h TU
       REAL*8 TSLORTO,TSLOCASO     ! TSL para orto y ocaso
       REAL*8 TUCREPM1,TUCREPM2   ! Tiempo universal de crepusculos matutino
       REAL*8 TUCREPV1,TUCREPV2   ! y vespertino utilizados para la iteracion
       REAL*8 TSLCREPM,TSLCREPV    ! TSL inicio crep. matutino y fin vespertino
       REAL XMIN,XMAX              ! Limites de la grafica en segundos.
       REAL YMIN,YMAX
       REAL YYY
       REAL XCREPV,XCREPM         ! Crepusculos en segundos.
       REAL*8 ARII(100),DECII(100)  ! AR y DEC del equinocio dado de entrada
       REAL*8 ARFF(100),DECFF(100)  ! AR y DEC de la fecha elegida
       REAL*8 AZI(100,100),ALT(100,100)! Azimut y altura  "     "
       REAL*8 TII,TFF               ! Equinocio de entrada y fecha (salida)
       REAL*8 MONTH
       REAL   TUPH(100),TUPS(100,100)!T.U. para plot estrella, en horas y seg
       REAL   TUPSU(100),ALTU(100)   ! T.U. y Alt para la ultima estrella.
       REAL   AZIU(100)              ! Azimut para la ultima estrella.
       REAL   SECZ(10)               ! Valores de sec Z, masa aire que utiliz.
       REAL   XAZI(100,100),YAZI(100,100) ! Valores de X e Y para el plot 
       REAL   RADIO(100)                  ! azimut-altura
       REAL   XMAZI(48),YMAZI(48)       ! Lo mismo para la luna
       REAL*8 NIGHTL        ! Duracion de la noche
       REAL*8 NIGHTLF       !   "       "  "  "     sin crepusculos
       REAL*8 ALTCUL,TUCUL  ! Altura y T.U. de la culminacion de un objeto
       REAL*8 TVIS(100)                 ! Tiempo que la estrella esta por 
       REAL*8 TVISINI(100),TVISFIN(100) ! encima horizonte en segundos 
       REAL*8 TUORSTAR,TUOCSTAR  ! T.U. del orto y el ocaso de una estrella. 
       REAL*8 DAYS,MALFA,MDELTA    
       REAL*8 MOONALT(48),MOONAZI(48) ! Altura y Azimut de Luna, cada hora   
       REAL*8 FASE(48)       ! Fase de la Luna
C
       CHARACTER*1  CH
       CHARACTER*2  CH2
       CHARACTER*11 CHDATE
       CHARACTER*13 LABEL(NWIN),VALUE(NWIN)
       CHARACTER*11 STARNAME(100),STARRA(100),STARDEC(100)
       CHARACTER*9  STARRAF(100),STARDECF(100)
       CHARACTER*6  EQUINOX(100)
       CHARACTER*10 OBSNAME(10),OBSLONG(10)
       CHARACTER*9  OBSLAT(10)     
       CHARACTER*5  OBSHEI(10)
       CHARACTER*11 TEXTBCK
       CHARACTER*1  ESTADO   ! Objeto circumpolar, invisible o con orto/ocaso  
       CHARACTER*7  CHORTOTU,CHOCASOTU    ! Hora en formato XXh XXm
       CHARACTER*7  CHORTOTSL,CHOCASOTSL  !  "       "         "
       CHARACTER*22 COPYR
       CHARACTER*3  CHSECZ(10)
       CHARACTER*11 INFILE     ! Listado de estrellas,N, Nombre, AR, DEC, Equin
       CHARACTER*41 LINE       ! con formato I2,A11,A9,A9,A6
       CHARACTER*6  CHORA
       CHARACTER*10 CFJ0
       CHARACTER*7  CNIGHTL,CNIGHTLF,CTFF
       CHARACTER*1  DIRCUL
       CHARACTER*5  CALTCUL
       CHARACTER*7  CTUCUL,CTVIS
       CHARACTER*4  CHF        ! Fase de la Luna
C     
       LOGICAL LELEV,LSECZ,LUT,LAZI
       LOGICAL LUNPLOT   ! TRUE ->  pgsci(0) -> borra el plot anterior
       LOGICAL LNSTAR    ! TRUE ->  Hay alguna estrella en el plot actual.
       LOGICAL LMOON(48) ! TRUE ->  La Luna es visible 
       LOGICAL VESIN     ! TRUE ->  Hay un fichero "vesinput" en ./ que se usara de input
C
       COMMON/BLKOBS/LAT,LONG
       COMMON/BLKOOC/HORTO,HOCASO,AORTO,AOCASO
c       COMMON/BLKDATE/FJ,IYEAR,IMONTH,IDAY,IHOUR,IMIN,ISEC
C
C----------------------------------------------------------------------
C     Definimos las ventanas
C
       DATA WIN /0.02,0.25,0.86,0.90,  ! 1  Star Name
     >           0.02,0.25,0.76,0.80,  ! 2  Ascension recta
     >           0.02,0.25,0.66,0.70,  ! 3  Declinacion
     >           0.02,0.13,0.56,0.60,  ! 4  Equinocio
     >           0.02,0.19,0.20,0.24,  ! 5  Draw
     >           0.02,0.19,0.10,0.14,  ! 6  Clear
     >           1.73,1.99,0.86,0.90,  ! 7  Nombre del Observatorio
     >           1.73,1.99,0.76,0.80,  ! 8  Longitud del Obs.
     >           1.73,1.99,0.66,0.70,  ! 9  Latitud del Obs.
     >           1.73,1.90,0.56,0.60,  ! 10 Altura del Obs.
     >           1.73,2.00,0.46,0.50,  ! 11 Fecha
     >           1.73,1.97,0.26,0.30,  ! 12 Crea listado estrellas, coord...
     >           1.73,1.97,0.18,0.22,  ! 13 Crea fichero postscript
     >           1.73,1.90,0.10,0.14,  ! 14 Exit
     >           0.30,0.36,0.50,0.60,  ! 15 Altura de la estrella
     >           0.30,0.36,0.30,0.40,  ! 16 Secante de Z, Masa de aire
     >           0.75,0.90,0.10,0.14,  ! 17 Tiempo Universal
     >           1.18,1.35,0.10,0.14,  ! 18 Azimut  
     >           0.60,1.50,0.00,0.06/  ! 19 I/O y ayuda.
C
       DATA LABEL /'Name','R.A.','Dec.','Equinox',' ',' ',
     >             'Observatory','Longitude','Latitude','Height','Date',
     >             ' ',' ',' ',' ',' ',' ',' ',' '/
C
       DATA VALUE /' -----',' 20 00 00','+20 00 00','2000',
     >             ' DRAW',' CLEAR',
     >             'Calar Alto','- 02 32 46','+37 13 25','2168m',' ',
     >             'file: ves.lis','file: ves.ps','  EXIT',
     >             'Elev.','sec z','  U.T.',' Azimut',' ' /
C
       DATA OBSNAME/'Calar Alto',' La Palma ',' La Silla ','  Madrid  ',
     >              ' Tenerife ','   Lick   ','  UKIRT   ','   KPNO   ',
     >              ' McDonald ','S.Pedro M.'/
C
       DATA OBSLONG/'- 02 32 46','- 17 52 47','- 70 44 00','- 03 41 15',
     >              '- 16 29 45','-123 38 12','-155 28 24','-111 36 00',
     >              '-104 01 18','-115 27 49'/
C
       DATA OBSLAT /'+37 13 25','+28 45 40','-29 15 00','+40 24 30',
     >              '+28 17 32','+37 20 36','+19 49 32','+31 57 48',
     >              '+30 40 18','+31 02 39'/
C
       DATA OBSHEI /'2168m','2334m','2400m',' 656m',
     >              '2400m','1283m','4194m','2120m',
     >              '2075m','2800m'/
C
       DATA SECZ /1.0,1.1,1.2,1.3,1.4,1.5,2.0,3.0,4.0,5.0/
       DATA CHSECZ /'1.0','1.1','1.2','1.3','1.4','1.5',
     >              '2.0','3.0','4.0','5.0'/
C
C----------------------------------------------------------------------
C
       NNOBS=1
       DO I=1,48
          LMOON(I)=.FALSE.
       END DO
C
       DEV(1)=PGOPEN('/XWINDOW')
       CALL PGSLCT(DEV(1))
       CALL PGPAGE       ! para que no se ponga negra la pantalla
       CALL PGSVP(0.0,1.0,0.0,1.0)
       CALL PGSWIN(0.0,2.0,0.0,1.0)
       CALL PGQWIN(X1,X2,Y1,Y2)
       CALL PGSCF(2)
       CALL PGSFS(1)
C
C----------------------------------------------------------------------
C     Le pedimos la fecha a la maquina
C
       CALL IDATE(INDATE)
       IDAY=INDATE(1)
       IMONTH=INDATE(2)
       IYEAR=INDATE(3)
       WRITE(CHDATE,"(I2,'-',I2,'-',I4)") IDAY,IMONTH,IYEAR
       WRITE(VALUE(11),'(A11)')CHDATE
       CALL ITIME(INTIME)
       IHOUR=INTIME(1)
       IMIN=INTIME(2)
       ISEC=INTIME(3)
       WRITE(CHORA,"(3(I2))")IHOUR,IMIN,ISEC
C
C----------------------------------------------------------------------
C     Dibujamos las ventanitas.
C
        DO I=1,NWIN
           CALL PGSCI(1)
           CALL PGRECT(WIN(1,I),WIN(2,I),WIN(3,I),WIN(4,I))
           X=WIN(1,I)
           Y=WIN(4,I)+0.01
           CALL PGSCH(0.8)
           CALL PGPTXT(X,Y,0.,0.0,LABEL(I))
           CALL PGSCI(4)
           CALL PGSCH(1.0)
           IF(I.EQ.15.OR.I.EQ.16) THEN
              ANGLE=90.
              Y=WIN(3,I)+0.02
              X=(WIN(1,I)+WIN(2,I)+0.02)/2.
           ELSE
              ANGLE=0.
              Y=WIN(3,I)+0.01
           END IF
           CALL PGPTXT(X,Y,ANGLE,0.0,VALUE(I))
        END DO 
C     La ventanita de I/O y ayuda tendra un fondo gris 
        CALL PGSCR(16,0.5,0.5,0.5)
        CALL IOWIN1('   Click center mouse button for HELP   ')
C
C     Al arrancar estaran activadas las opciones de Elevacion en eje Y
C     y Tiempo Universal en eje X.
        CALL PGSCI(4)
        CALL PGSFS(2)
        CALL PGSLW(12)
        CALL PGRECT(WIN(1,15),WIN(2,15),WIN(3,15),WIN(4,15))
        CALL PGRECT(WIN(1,17),WIN(2,17),WIN(3,17),WIN(4,17))
        CALL PGSLW(1)
        CALL PGQWIN(X1,X2,Y1,Y2)
        LELEV=.TRUE.
        LSECZ=.FALSE.
        LUT=.TRUE.
        LAZI=.FALSE.
        LNSTAR=.FALSE.    ! Al empezar no hay ninguna estrella en el plot
C
C     El observatorio inicial es Calar Alto, pasamos sus coordenadas
C     a valor decimal.
        CALL CH2DEC(VALUE(8),LONG,1)
        LONG=LONG/15.D0              ! Pasamos de grados a horas.
        CALL CH2DEC(VALUE(9),LAT,2)
        READ(VALUE(10)(1:4),*) ALTOBS
C
        WRITE(COPYR,'(A22)')'\(0274) pedraz@caha.es'
        CALL PGSCI(1)
        CALL PGSCH(0.55)
        CALL PGPTXT(1.99,0.01,0.,1.,COPYR)
        CALL PGSCH(1.0)
C
C----------------------------------------------------------------------
C     Calculamos los momentos de anochecer, amanecer y crepusculos para
C     delimitar la grafica.
C 
 50     LUNPLOT=.FALSE.
C
C- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
C     Calculamos la Fecha Juliana (FJ) y TS0 a las 0h TU en Greenwich. 
C     A partir de estas se calcula la posicion del sol. Iteramos para 
C     obtener la posicion del sol en el ocaso. Como valor inicial 
C     tomamos la posicion a las 0h TU
      CALL FJ_TS0(IYEAR,IMONTH,IDAY,0.D0,FJ,TS0)
       FJ0=FJ
       TS1=TS0+24.*2.737909E-3
       IF(TS1.GT.24.D0) TS1=TS1-24.D0
       CALL POSSOL(FJ,ARSUN,DECSUN)
C
C     Calculamos la distancia zenital en el ocaso teniendo en cuenta la
C     refraccion (34'), semidiametro del Sol (16') y la depresion del
C     horizonte(elipsoide WGS-84) dada la altura del observatorio.
C
       SEMIA=6378.137D0
       SEMIB=6356.752D0
       REARTH=1.D0/
     >        SQRT(COS(LAT*PI/180.D0)*COS(LAT*PI/180.D0)/SEMIA/SEMIA+
     >        SIN(LAT*PI/180.D0)*SIN(LAT*PI/180.D0)/SEMIB/SEMIB)
       ANGHOA=ACOS(REARTH/(REARTH+ALTOBS/1000.D0))
       ANGHOA=ANGHOA*180.D0/PI
       DISZEN=90.D0+34.D0/60.D0+16.D0/60.D0+ANGHOA
C
C- - - - - - - - - - -  OCASO - - - - - -  - - - - - - - - - - - - - 
C     Con la declinacion y la distancia zenital calculamos angulo
C     horario y azimut del sol para el ocaso
      CALL ORCASO(DECSUN,DISZEN,ESTADO)
      IF(ESTADO.EQ.'C')THEN
        CALL IOWIN1(' Circumpolar SUN  --  NO orto-ocaso     ')
        GOTO 70
      ELSE IF(ESTADO.EQ.'I')THEN
	CALL IOWIN1(' Invisible SUN  --  NO orto-ocaso       ')
        GOTO 70
      END IF
C     Iterando, calculamos el ocaso que se produce despues de las 0h TU
C     del dia de la observacion.
      CALL TIEMPOU(ARSUN,HOCASO,LONG,TS0,TSL,TUOCASO1)
 52   CALL FJ_TS0(IYEAR,IMONTH,IDAY,TUOCASO1,FJ,TS0)
      CALL POSSOL(FJ,ARSUN,DECSUN)
      CALL ORCASO(DECSUN,DISZEN,ESTADO)
      CALL TIEMPOU(ARSUN,HOCASO,LONG,TS0,TSL,TUOCASO2)
      TSLOCASO=TSL
      IF(ABS(TUOCASO1-TUOCASO2).GT.0.0003D0) THEN
         TUOCASO1=TUOCASO2
         GOTO 52
      END IF
C
C - - - - - - -  Crepusculo ASTRONOMICO vespertino - - - - - - - - -
C     Como valor inicial tomamos el ocaso.
       TUCREPV1=TUOCASO2
C
 53    CALL FJ_TS0(IYEAR,IMONTH,IDAY,TUCREPV1,FJ,TS0)
       CALL POSSOL(FJ,ARSUN,DECSUN)
       CALL ORCASO(DECSUN,108.0+ANGHOA,ESTADO)
       IF(ESTADO.EQ.'C') THEN
          CALL IOWIN1(' Twilight all night long                ')
          GOTO 70
       END IF
       CALL TIEMPOU(ARSUN,HOCASO,LONG,TS0,TSL,TUCREPV2)
       TSLCREPV=TSL
       IF(TUCREPV2.LT.TUCREPV1) TUCREPV2=TUCREPV2+24.0
       IF(ABS(TUCREPV2-TUCREPV1).GT.0.0003D0) THEN
         TUCREPV1=TUCREPV2
         GOTO 53
       END IF
C
C - - - - - - - - - -  ORTO - - - - - - - - - - - - - - - - - - - - 
C     Calculamos el orto del dia siguiente. Partimos de que a las 0h TU
C     la fecha juliana sera un dia mas FJ1=FJ0+1.
       FJ1=FJ0+1.D0
       IDAY=IDAY+1
 54    CALL POSSOL(FJ1,ARSUN,DECSUN)
       CALL ORCASO(DECSUN,DISZEN,ESTADO)
       CALL TIEMPOU(ARSUN,HORTO,LONG,TS0,TSL,TUORTO1)
C     Iteramos
 55    CALL FJ_TS0(IYEAR,IMONTH,IDAY,TUORTO1,FJ,TS0)
       CALL POSSOL(FJ,ARSUN,DECSUN)
       CALL ORCASO(DECSUN,DISZEN,ESTADO)
       CALL TIEMPOU(ARSUN,HORTO,LONG,TS0,TSL,TUORTO2)
       TSLORTO=TSL
       IF(ABS(TUORTO1-TUORTO2).GT.0.0003D0) THEN
          TUORTO1=TUORTO2
          GOTO 55
       END IF
C
C - - - - - - -  Crepusculo ASTRONOMICO matutino - - - - - - - - -
C     Como valor inicial tomamos el orto.
      TUCREPM1=TUORTO2
C
 56    CALL FJ_TS0(IYEAR,IMONTH,IDAY,TUCREPM1,FJ,TS0)
       CALL POSSOL(FJ,ARSUN,DECSUN)
       CALL ORCASO(DECSUN,108.0+ANGHOA,ESTADO)
       CALL TIEMPOU(ARSUN,HORTO,LONG,TS0,TSL,TUCREPM2)
       TSLCREPM=TSL
       IF(TUCREPM2.GT.TUCREPM1) TUCREPM2=TUCREPM2-24.0
       IF(ABS(TUCREPM2-TUCREPM1).GT.0.0003D0) THEN
          TUCREPM1=TUCREPM2
          GOTO 56
       END IF
       IF(FJ1.EQ.FJ0) GOTO 57
       IDAY=IDAY-1
       IF(TUORTO2.GT.TUOCASO2) THEN ! Puede ser que debido a su longitud
          FJ1=FJ0                   ! en un observatorio anochezca y 
          GOTO 54                   ! amanezca en el mismo dia juliano.
       END IF
 57    CONTINUE
C
C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  
C     Los limites de la grafica seran el comienzo de la noche
C     ocaso=TUOCASO2 y el final orto=TUORTO2. El final del crepusculo
C     por la tarde es TUCREPV2 y el inicio del crepusculo por la
C     maniana TUCREPM2.
C
       XMIN=TUOCASO2*3600.
       XCREPV=TUCREPV2*3600.
       IF(XCREPV.LT.XMIN) XCREPV=XCREPV+(24.*3600.)
       IF(FJ1.EQ.FJ0) THEN 
          XMAX=TUORTO2*3600.
          XCREPM=TUCREPM2*3600.
       ELSE
          XMAX=(TUORTO2+24.)*3600.
          XCREPM=(TUCREPM2+24.)*3600.
       END IF 
       YMIN=0.
       YMAX=90.
 60    CALL PGSVP(0.22,0.85,0.22,0.90) !Definimos el area del plot que 
C                                      !utilizaremos para la grafica   
       CALL PGSWIN(XMIN,XMAX,YMIN,YMAX)
C
       IF(LUT) THEN
          IF(LELEV) THEN
             CALL PGSCI(0)         ! Borramos
             CALL PGSLS(1) 
             CALL PGTBOX('ABCNTSZYXHO',0.0,0,'ABC',0.0,0)
             CALL PGSLS(4) 
             DO I=1,10              
                Y=ACOS(1./SECZ(I)) ! Lineas horizontales para marcar las
                Y=Y*180./PI        ! sec Z (masas de aire)
                Y=90.-Y
                CALL PGMOVE(XMIN,Y)
                CALL PGDRAW(XMAX,Y)
                CALL PGPTXT(XMIN-300.,Y-1.,0.0,1.0,CHSECZ(I))
             END DO
             CALL PGSCI(1)      ! Pintamos
             IF(LUNPLOT) CALL PGSCI(0)
             CALL PGSLS(1) 
             CALL PGTBOX('ABCNTSZYXHO',0.0,0,'ABCNTS',0.0,0)
             CALL PGSLS(4) 
             DO I=10,90,10      ! lineas horizontales cada 10 grados de elevac.
                CALL PGMOVE(XMIN,REAL(I))
                CALL PGDRAW(XMAX,REAL(I))
             END DO
          ELSE IF(LSECZ) THEN
             CALL PGSCI(0)      ! Borramos
             CALL PGSLS(1) 
             CALL PGTBOX('ABCNTSZYXHO',0.0,0,'ABCNTS',0.0,0)
             CALL PGSLS(4) 
             DO I=10,90,10      ! lineas horizontales cada 10 grados de elevac.
                CALL PGMOVE(XMIN,REAL(I))
                CALL PGDRAW(XMAX,REAL(I))
             END DO
             CALL PGSCI(1)         ! Pintamos
             IF(LUNPLOT) CALL PGSCI(0)
             CALL PGSLS(1) 
             CALL PGTBOX('ABCNTSZYXHO',0.0,0,'ABC',0.0,0)
             CALL PGSLS(4) 
             DO I=1,10              
                Y=ACOS(1./SECZ(I)) ! Lineas horizontales para marcar las
                Y=Y*180.D0/PI        ! sec Z (masas de aire)
                Y=90.-Y
                CALL PGMOVE(XMIN,Y)
                CALL PGDRAW(XMAX,Y)
                CALL PGPTXT(XMIN-300.,Y-1.,0.0,1.0,CHSECZ(I))
             END DO
          END IF
          IF(LUNPLOT) CALL PGSCI(0)
          CALL PGSLS(4)
          DO I=0,172800,3600        ! Lineas de puntos verticales cada hora
             IF(I.GT.XMIN.AND.I.LT.XMAX) THEN
                CALL PGMOVE(REAL(I),0.)
                CALL PGDRAW(REAL(I),90.)
             END IF
          END DO
          CALL PGSLS(2)  ! Lineas verticales a trazos para marcar los crepusc
          CALL PGMOVE(XCREPV,0.)
          CALL PGDRAW(XCREPV,90.)
          CALL PGMOVE(XCREPM,0.)
          CALL PGDRAW(XCREPM,90.)
          CALL PGSLS(1)
       END IF
C     Labels con tiempos de orto y ocaso en TU y TSL
       IF(LUNPLOT) THEN
          CALL PGSCI(0)
       ELSE
          CALL PGSCI(1)
       END IF
       CALL PGSCH(0.8)
       CALL PGPTXT(XMIN,91.,0.0,0.0,'Sunset')
       CALL R2CH(TUOCASO2,CHOCASOTU)
       CALL PGPTXT(XMIN,94.,0.0,0.0,CHOCASOTU//' UT')
       CALL R2CH(TSLOCASO,CHOCASOTSL)
       CALL PGPTXT(XMIN,97.,0.0,0.0,CHOCASOTSL//' LST')
       CALL PGPTXT(XMAX,91.,0.0,1.0,'Sunrise')
       CALL R2CH(TUORTO2,CHORTOTU)
       CALL PGPTXT(XMAX,94.,0.0,1.0,CHORTOTU//'  UT ')
       CALL R2CH(TSLORTO,CHORTOTSL)
       CALL PGPTXT(XMAX,97.,0.0,1.0,CHORTOTSL//' LST')
       CALL PGSCH(1.0)
C
       IF(LUNPLOT) GOTO 50
C
C----------------------------------------------------------------------
C     Calculamos la posicion de la Luna para cada hora 
C
       TSL0=TS0+LONG   ! Tiempo sidereo local a las 0h TU
C
       DO I=1,48
          CALL NDAYS(IYEAR,IMONTH,IDAY,I,0,0,DAYS) !Dias transcurr del anio 
          CALL MOONPOS(FJ0,IYEAR,DAYS,MALFA,MDELTA,FASE(I))
          TSL=TSL0+DBLE(I)+DBLE(I)*2.737909D-3
          IF(TSL.GT.24.D0) TSL=TSL-24.D0
          IF(TSL.GT.24.D0) TSL=TSL-24.D0
          CALL ECU2HOR(TSL,LAT,MALFA,MDELTA,MOONAZI(I),MOONALT(I))
          IF(MOONALT(I).GT.0.D0) THEN
             IF(I*3600.GT.XMIN.AND.I*3600.LT.XMAX) THEN
                LMOON(I)=.TRUE.
                RADIO(I)=(90.-MOONALT(I))/100.D0
                IF(MOONAZI(I).GT.90.D0) THEN
                   MOONAZI(I)=360.D0+90.D0-MOONAZI(I)
                ELSE
                   MOONAZI(I)=90.D0-MOONAZI(I)
                END IF
                XMAZI(I)=RADIO(I)*COS(MOONAZI(I)*PI/180.D0)
                YMAZI(I)=RADIO(I)*SIN(MOONAZI(I)*PI/180.D0)
             ELSE
                LMOON(I)=.FALSE.
             END IF
          ELSE
             LMOON(I)=.FALSE.
          END IF
       END DO
C
       CALL PGSCI(5)
       IF(LUT) THEN
          DO I=1,48
             IF(LMOON(I)) THEN 
                CALL PGSCH(1.5)
                YYY=MOONALT(I)
                CALL PGPT1(I*3600.,YYY,2291)
                WRITE(CHF,'(F4.2)')FASE(I)
                CALL PGSCH(0.6)
                CALL PGTEXT(I*3600.,YYY-0.5,CHF)
             END IF
          END DO
       ELSE IF(LAZI) THEN
          CALL PGSVP(0.24,0.82,0.20,0.94)
          CALL PGSWIN(-1.,1.,-1.,1.)
          DO I=1,48
             IF(LMOON(I)) THEN
                CALL PGSCH(1.5)
                CALL PGPT1(XMAZI(I),YMAZI(I),2291)
                WRITE(CHF,'(F4.2)')FASE(I)
                CALL PGSCH(0.6)
                CALL PGTEXT(XMAZI(I),YMAZI(I)-0.01,CHF)
             END IF
          END DO
       END IF
C
       CALL PGSCI(1)
       CALL PGSCH(1.)
C
C----------------------------------------------------------------------
C      Miramos en el directorio ./ a ver si hay un fichero "vesinput"
C      que se utilizara como lista de objetos para pintar en la 
C      grafica al lanzar ves
       INQUIRE(FILE='vesinput',EXIST=VESIN)
C----------------------------------------------------------------------
C     Loop para elegir ventanita con el cursor.
C
       IF(LNSTAR) GOTO 500
C
        XX=1.0
        YY=0.5
        NSTAR=0
C
 70     IF(VESIN) GOTO 100
        CALL PGSVP(0.0,1.0,0.0,1.0)
        CALL PGSWIN(0.0,2.0,0.0,1.0)
        CALL PGSLW(8)
        CALL PGSCI(4)
        DO 90 N=1,1000
           CALL PGBAND(0,0,XX,YY,XX,YY,CH)           
           IF(ICHAR(CH).EQ.0) GOTO 2000
           CALL PGQWIN(X1,X2,Y1,Y2)
           DO 80 NN=1,NWIN
              IF (WIN(1,NN).LE.XX .AND. WIN(2,NN).GE.XX .AND.
     >            WIN(3,NN).LE.YY .AND. WIN(4,NN).GE.YY) THEN
                 IF(CH.EQ.'D') THEN
                    CALL HELP(NN)
                    GOTO 70
                 END IF
                 IF(NN.EQ.1) GOTO 100
                 IF(NN.EQ.2) GOTO 200
                 IF(NN.EQ.3) GOTO 300
                 IF(NN.EQ.4) GOTO 400
                 IF(NN.EQ.5) GOTO 500
                 IF(NN.EQ.6) THEN
                    LNSTAR=.FALSE.
                    GOTO 600
                 END IF
                 IF(NN.EQ.7) GOTO 700
                 IF(NN.EQ.8) GOTO 800
                 IF(NN.EQ.9) GOTO 900
                 IF(NN.EQ.10) GOTO 1000
                 IF(NN.EQ.11) GOTO 1100
                 IF(NN.EQ.12) GOTO 1200
                 IF(NN.EQ.13) GOTO 1300
                 IF(NN.EQ.14) GOTO 1400 
                 IF(NN.EQ.15) GOTO 1500 
                 IF(NN.EQ.16) GOTO 1600 
                 IF(NN.EQ.17) GOTO 600 
                 IF(NN.EQ.18) GOTO 600
              END IF
 80     CONTINUE
              CALL IOWIN1('No window selected                      ')
 90    CONTINUE   
C----------------------------------------------------------------------
C     Si pinchamos con el boton izquierdo introducimos el nombre de una 
C     nueva estrella. Si es el boton derecho podemos dar el nombre de 
C     un fichero que contenga varias estrellas.
C
 100   IF(VESIN) WRITE(CH,'(A1)')'X'
       IF(CH.EQ.'A') THEN         ! Boton izquierdo del raton
          CALL PGRECT(WIN(1,1),WIN(2,1),WIN(3,1),WIN(4,1))
          CALL PGRSTR(WIN(1,1),WIN(3,1)+0.01,VALUE(1),TEXTBCK,NN,0)
          CALL PGSCI(1)
          CALL PGSLW(8)
          CALL PGRECT(WIN(1,1),WIN(2,1),WIN(3,1),WIN(4,1))
          CALL PGSLW(8)
          GOTO 200
       ELSE IF(CH.EQ.'D') THEN    ! Boton central del raton
          CALL IOWIN1('AYUDA                                   ')
       ELSE IF(CH.EQ.'X') THEN    ! Boton derecho del raton
          IF(VESIN) THEN
            WRITE(INFILE,'(A11)') 'vesinput   '
            VESIN=.FALSE.
            GOTO 160
          END IF
       CALL IOWIN2('File format :                                     ',
     >             'Name(4:14),R.A.(16:24),DEC.(26:34),EQUI.(36:41)   ')
          GOTO 150
 140      CALL IOWIN1('No file with this name, type again ...  ')
 150      CALL PGSCI(4)
          CALL PGSFS(2)
          CALL PGSLW(8)
          CALL PGRECT(WIN(1,1),WIN(2,1),WIN(3,1),WIN(4,1))
          CALL PGRSTR(WIN(1,1),WIN(3,1)+0.01,VALUE(1),TEXTBCK,NN,0)
          CALL PGSCI(1)
          CALL PGSLW(8)
          CALL PGRECT(WIN(1,1),WIN(2,1),WIN(3,1),WIN(4,1))
          CALL PGSLW(1)
          VALUE(1)=TEXTBCK
          INFILE=TEXTBCK
          IF(INFILE.EQ.' ') THEN
             CALL IOWIN1(' No filename typed  ... OK              ')
             GOTO 70
          END IF
 160      OPEN(21,FILE=INFILE,STATUS='OLD',FORM='FORMATTED',ERR=140)
          DO I=1,100
             NSTAR=NSTAR+1
             READ(21,'(A41)',END=190) LINE
             STARNAME(NSTAR)=LINE(4:14)
             STARRA(NSTAR)=LINE(16:24)
             STARDEC(NSTAR)=LINE(26:34)
             EQUINOX(NSTAR)=LINE(36:41)
             CALL CH2DEC(STARRA(I),ARII(NSTAR),2)
             CALL CH2DEC(STARDEC(I),DECII(NSTAR),2)
          END DO
          CALL IOWIN1('    More than 100 stars in this file    ') 
       ELSE
          GOTO 70
       END IF
 190   CLOSE(21)
       NSTAROLD=NSTAR-1
       NSTAR=0
       GOTO 500
C----------------------------------------------------------------------
C     Rellenamos el valor de la ascension recta de la estrella.
C
 200   CALL PGSCI(4)
       CALL PGRECT(WIN(1,2),WIN(2,2),WIN(3,2),WIN(4,2))
       CALL PGRSTR(WIN(1,2),WIN(3,2)+0.01,VALUE(2),TEXTBCK,NN,1)
       VALUE(2)=TEXTBCK
       CALL PGSCI(1)
       CALL PGSLW(8)
       CALL PGRECT(WIN(1,2),WIN(2,2),WIN(3,2),WIN(4,2))
C
       IF(VALUE(2)(1:1).NE.'-'.AND.VALUE(2)(1:1).NE.'+'
     >                        .AND.VALUE(2)(1:1).NE.' ') GOTO 250
       IF(VALUE(2)(4:4).NE.' '.OR.VALUE(2)(7:7).NE.' ') GOTO 250
       DO J=2,9
          IF(TESTN(VALUE(2)(J:J)).EQ.0.AND.
     >       VALUE(2)(J:J).NE.' ') GOTO 250 
       END DO
C
       IF(NN.EQ.1) THEN 
          GOTO 300
       ELSE
          GOTO 70
       END IF
C
 250   CALL IOWIN1(' Please use format :  00 00 00          ')
       GOTO 200
C----------------------------------------------------------------------
C     Rellenamos el valor de la declinacion de la estrella.
C
 300   CALL PGSCI(4)
       CALL PGRECT(WIN(1,3),WIN(2,3),WIN(3,3),WIN(4,3))
       CALL PGRSTR(WIN(1,3),WIN(3,3)+0.01,VALUE(3),TEXTBCK,NN,0)
       VALUE(3)=TEXTBCK
       CALL PGSCI(1)
       CALL PGSLW(8)
       CALL PGRECT(WIN(1,3),WIN(2,3),WIN(3,3),WIN(4,3))
C
       IF(VALUE(3)(1:1).NE.'-'.AND.VALUE(3)(1:1).NE.'+'
     >                        .AND.VALUE(3)(1:1).NE.' ') THEN
C         ! Borramos el texto
          CALL PGPTXT (WIN(1,3),WIN(3,3)+0.01,0.0,0.0,VALUE(3))
          DO J=9,2,-1
             VALUE(3)(J:J)=VALUE(3)(J-1:J-1)
          END DO
          VALUE(3)(1:1)='+'
          CALL PGSCI(4)
          CALL PGSLW(1)
          CALL PGPTXT (WIN(1,3),WIN(3,3)+0.01,0.0,0.0,VALUE(3))
          CALL PGSCI(1)
       END IF
       IF(VALUE(3)(4:4).NE.' '.OR.VALUE(3)(7:7).NE.' ') GOTO 350
       DO J=2,9
          IF(TESTN(VALUE(3)(J:J)).EQ.0.AND.
     >       VALUE(3)(J:J).NE.' ') GOTO 350 
       END DO
C
       IF(NN.EQ.1) THEN 
          NN=4
          GOTO 400
       ELSE
          GOTO 70
       END IF
C
 350   CALL IOWIN1('Please use format :  +00 00 00          ')
       GOTO 300
C----------------------------------------------------------------------
C     Rellenamos el valor del equinocio para la A.R y DEC de la nueva estrella
C
 400   CALL PGSCI(4)
       CALL PGRECT(WIN(1,4),WIN(2,4),WIN(3,4),WIN(4,4))
       CALL PGRSTR(WIN(1,4),WIN(3,4)+0.01,VALUE(4),TEXTBCK,NN,0)
       CALL PGSCI(1)
       CALL PGSLW(8)
       CALL PGRECT(WIN(1,4),WIN(2,4),WIN(3,4),WIN(4,4))
       CALL PGSLW(1)
       VALUE(4)=TEXTBCK
C
       DO J=1,4
          IF(TESTN(VALUE(4)(J:J)).EQ.0.AND.
     >       VALUE(4)(J:J).NE.' ') GOTO 450 
       END DO
C
       GOTO 70
C
 450   CALL IOWIN1('Please use format :  0000               ')
       GOTO 400
C----------------------------------------------------------------------
C     Trazamos la posicion de la nueva estrella que tiene las coordenadas
C     que se ven en las correspondientes ventanas. La grafica corresponde
C     a los ejes que estan resaltados.
 500   LNSTAR=.TRUE.
       IF(NN.EQ.5) THEN
          CALL PGSLW(12)
          CALL PGSCI(4)
          CALL PGRECT(WIN(1,5),WIN(2,5),WIN(3,5),WIN(4,5))
          NSTAR=NSTAR+1
          STARNAME(NSTAR)=VALUE(1)
          STARRA(NSTAR)=VALUE(2)
          STARDEC(NSTAR)=VALUE(3)
          EQUINOX(NSTAR)=VALUE(4)
          CALL CH2DEC(VALUE(2),ARII(NSTAR),2)
          CALL CH2DEC(VALUE(3),DECII(NSTAR),2)
       END IF
C
       READ(VALUE(4)(1:4),*) TII
       IF(IMONTH.LE.2)THEN
         MONTH=AINT((IMONTH-1)*63./2.)
       ELSE
         MONTH=AINT((IMONTH+1)*30.6)-63.
       END IF
       MONTH=MONTH+IDAY
       TFF=REAL(IYEAR+MONTH/365.25)
       WRITE(CTFF,'(F7.2)') TFF
C
       IF(NSTAR.EQ.0) THEN ! Hemos cambiado alguna ventana de la 7 a la 11
          CALL PGSCI(4)    ! o hemos metido un listado de estrellas nuevas
          CALL PGSLW(2)
          DO J=1,NSTAROLD
c             TVISINI(J)=DREAL(XMIN)
c             TVISFIN(J)=DREAL(XMAX)
             TVISINI(J)=XMIN
             TVISFIN(J)=XMAX
             CALL CHANEQUI(TII,TFF,ARII(J),ARFF(J),
     >                     DECII(J),DECFF(J))
             DO I=1,100
                TUPS(J,I)=XMIN+(XMAX-XMIN)*REAL(I-1)/99.
                TUPH(I)=TUPS(J,I)/3600.
                TSL=TS0+TUPH(I)+TUPH(I)*2.737909E-3+LONG
                CALL ECU2HOR(TSL,LAT,ARFF(J),DECFF(J),
     >                       AZI(J,I),ALT(J,I))
                TUPSU(I)=TUPS(J,I)
                ALTU(I)=ALT(J,I)
                AZIU(I)=AZI(J,I)
                RADIO(I)=(90.-ALTU(I))/100.
                IF(AZIU(I).GT.90.) THEN
                   AZIU(I)=360.+90.-AZIU(I)
                ELSE
                   AZIU(I)=90.-AZIU(I)
                END IF
                XAZI(J,I)=RADIO(I)*COS(AZIU(I)*PI/180.)
                YAZI(J,I)=RADIO(I)*SIN(AZIU(I)*PI/180.)
             END DO
             WRITE(CH2,'(I2)') J
             CALL PGSCH(0.5)
             IF(LUT) THEN 
                CALL PGSVP(0.22,0.85,0.22,0.90)
                CALL PGSWIN(XMIN,XMAX,YMIN,YMAX)
                CALL PGLINE(100,TUPSU,ALTU)
                DO I=1,9
                   IF(ALTU(10*I+J).GT.0.AND.ALTU(10*I+J).LT.89) THEN
                      CALL PGTEXT(TUPSU(10*I+J),ALTU(10*I+J),CH2)
                   END IF
                END DO
             ELSE IF(LAZI) THEN
                CALL PGSVP(0.24,0.82,0.20,0.94)
                CALL PGSWIN(-1.,1.,-1.,1.)
                DO I=1,100
                   IF(I.GT.1.AND.ALTU(I-1).GT.0.AND.ALTU(I).GT.0.)THEN
                      CALL PGMOVE(XAZI(J,I-1),YAZI(J,I-1))
                      CALL PGDRAW(XAZI(J,I),YAZI(J,I))
                   END IF
                END DO
                DO I=1,9
                   IF(ALTU(10*I+J).GT.0.AND.ALTU(10*I+J).LT.89) THEN
                      CALL PGTEXT(XAZI(J,10*I+J),YAZI(J,10*I+J),CH2)
                   END IF
                END DO
             END IF
             CALL PGSCH(1.0)
          END DO
          CALL PGSVP(0.0,1.0,0.0,1.0)
          CALL PGSWIN(0.0,2.0,0.0,1.0)
          NSTAR=NSTAROLD
          GOTO 70
       ELSE
c          TVISINI(NSTAR)=DREAL(XMIN)
c          TVISFIN(NSTAR)=DREAL(XMAX)
          TVISINI(NSTAR)=XMIN
          TVISFIN(NSTAR)=XMAX
          DO I=1,100
             TUPS(NSTAR,I)=XMIN+(XMAX-XMIN)*REAL(I-1)/99.
             TUPH(I)=TUPS(NSTAR,I)/3600.
             TSL=TS0+TUPH(I)+TUPH(I)*2.737909D-3+LONG
             CALL CHANEQUI(TII,TFF,ARII(NSTAR),ARFF(NSTAR),
     >                     DECII(NSTAR),DECFF(NSTAR))
             CALL ECU2HOR(TSL,LAT,ARFF(NSTAR),DECFF(NSTAR),
     >                    AZI(NSTAR,I),ALT(NSTAR,I))
             TUPSU(I)=TUPS(NSTAR,I)
             ALTU(I)=ALT(NSTAR,I)
             AZIU(I)=AZI(NSTAR,I)
             IF(ALTU(I).GT.0.) THEN
                RADIO(I)=(90.-ALTU(I))/100.
                IF(AZIU(I).GT.90.) THEN
                   AZIU(I)=360.+90.-AZIU(I)
                ELSE
                   AZIU(I)=90.-AZIU(I)
                END IF
                XAZI(NSTAR,I)=RADIO(I)*COS(AZIU(I)*PI/180.)
                YAZI(NSTAR,I)=RADIO(I)*SIN(AZIU(I)*PI/180.)
             END IF
          END DO
C
          CALL PGSCI(4)
          CALL PGSLW(2)
          CALL PGSVP(0.22,0.85,0.22,0.90)
          WRITE(CH2,'(I2)') NSTAR
          CALL PGSCH(0.5)
          IF(LUT) THEN
             CALL PGSVP(0.22,0.85,0.22,0.90)
             CALL PGSWIN(XMIN,XMAX,YMIN,YMAX)
             CALL PGLINE(100,TUPSU,ALTU)
             DO I=1,8
                IF(ALTU(10*I+NSTAR).GT.0.) THEN
                   CALL PGTEXT(TUPSU(10*I+NSTAR),ALTU(10*I+NSTAR),CH2)
                END IF
             END DO
          ELSE IF(LAZI) THEN
             CALL PGSVP(0.24,0.82,0.20,0.94)
             CALL PGSWIN(-1.,1.,-1.,1.)
             DO I=1,100
                   IF(I.GT.1.AND.ALTU(I-1).GT.0.AND.ALTU(I).GT.0.)THEN
                      CALL PGMOVE(XAZI(NSTAR,I-1),YAZI(NSTAR,I-1))
                      CALL PGDRAW(XAZI(NSTAR,I),YAZI(NSTAR,I))
                   END IF
             END DO
             DO I=1,9
                IF(ALTU(10*I+J).GT.0.) THEN
                  CALL PGTEXT(XAZI(NSTAR,10*I+J),YAZI(NSTAR,10*I+J),CH2)
                END IF
             END DO
          END IF
          CALL PGSCH(1.0)
          CALL PGSVP(0.0,1.0,0.0,1.0)
          CALL PGSWIN(0.0,2.0,0.0,1.0)
       END IF
C
       CALL PGSCI(1)
       CALL PGSLW(12)
       CALL PGRECT(WIN(1,5),WIN(2,5),WIN(3,5),WIN(4,5))
       CALL PGSLW(1)
       GOTO 70

C----------------------------------------------------------------------
C     Quitamos las trazas de las estrellas calculadas hasta ahora.
C
 600   CALL PGSLW(12)
       CALL PGSCI(4)
       CALL PGSLS(1)
       CALL PGRECT(WIN(1,6),WIN(2,6),WIN(3,6),WIN(4,6))
C
       CALL PGSCI(0)
       CALL PGSLW(2)
       CALL PGSVP(0.22,0.85,0.22,0.90)
       CALL PGSWIN(XMIN,XMAX,YMIN,YMAX)
       DO J=1,NSTAR
          DO I=1,100
             TUPSU(I)=TUPS(J,I)
             ALTU(I)=ALT(J,I)
          END DO
          WRITE(CH2,'(I2)') J
          CALL PGSCH(0.5)
          IF(LUT) THEN
             CALL PGLINE(100,TUPSU,ALTU)
             DO I=1,9
                IF(ALTU(10*I+J).GT.0.AND.ALTU(10*I+J).LT.89) THEN
                   CALL PGTEXT(TUPSU(10*I+J),ALTU(10*I+J),CH2)
                END IF
             END DO
          ELSE 
             CALL PGSVP(0.24,0.82,0.20,0.94)
             CALL PGSWIN(-1.,1.,-1.,1.)
             DO I=1,100
                CALL PGDRAW(XAZI(J,I),YAZI(J,I))
             END DO
             DO I=1,9
                IF(ALTU(10*I+J).GT.0.AND.ALTU(10*I+J).LT.89) THEN
                   CALL PGTEXT(XAZI(J,10*I+J),YAZI(J,10*I+J),CH2)
                END IF
             END DO
          END IF
          CALL PGSCH(1.0)
       END DO
C     Borramos la Luna
       IF(LUT) THEN
          DO I=1,48
             IF(LMOON(I)) THEN
                CALL PGSCH(1.5)
                YYY=MOONALT(I)
                CALL PGPT1(I*3600.,YYY,2291)
                WRITE(CHF,'(F4.2)')FASE(I)
                CALL PGSCH(0.6)
                CALL PGTEXT(I*3600.,YYY-0.5,CHF)
             END IF
          END DO
       ELSE IF(LAZI) THEN
          CALL PGSVP(0.24,0.82,0.20,0.94)
          CALL PGSWIN(-1.,1.,-1.,1.)
          DO I=1,48
             IF(LMOON(I)) THEN
                CALL PGSCH(1.5)
                CALL PGPT1(XMAZI(I),YMAZI(I),2291)
                WRITE(CHF,'(F4.2)')FASE(I)
                CALL PGSCH(0.6)
                CALL PGTEXT(XMAZI(I),YMAZI(I)-0.01,CHF)
             END IF
          END DO
       END IF
       CALL PGSCH(1.)
       IF(NN.EQ.6) THEN
          DO I=1,48
             LMOON(I)=.FALSE.
          END DO
       END IF
C
       CALL PGSLW(12)
       CALL PGSCI(1)
       CALL PGSVP(0.0,1.0,0.0,1.0)
       CALL PGSWIN(0.0,2.0,0.0,1.0)
       CALL PGRECT(WIN(1,6),WIN(2,6),WIN(3,6),WIN(4,6))
       CALL PGSLW(1)

C
       IF(NN.EQ.17) GOTO 1700
       IF(NN.EQ.18) GOTO 1800
       NSTAROLD=NSTAR
       NSTAR=0
       IF(NN.EQ.7.AND.CH.EQ.'A') GOTO 800
       IF(NN.EQ.7.AND.CH.EQ.'X') GOTO 60
       IF(NN.GE.8.AND.NN.LE.11) GOTO 60
       GOTO 70
C
C----------------------------------------------------------------------
C     Si pinchamos con el boton izquierdo podemos meter los datos de un
C     nuevo observatorio. Hemos de rellenar las cuatro ventanas corres-
C     pondientes. Nombre, longitud, latitud y altura.
C     Si pinchamos con el boton derecho vamos elegiendo otro observatorio
C     de los diez que hay en el programa (DATA).
C
 700   IF(CH.EQ.'A') THEN         ! Boton izquierdo del raton
          CALL PGRECT(WIN(1,7),WIN(2,7),WIN(3,7),WIN(4,7))
          CALL PGRSTR(WIN(1,7),WIN(3,7)+0.01,VALUE(7),TEXTBCK,NN,0)
          CALL PGSCI(1)
          CALL PGSLW(8)
          CALL PGRECT(WIN(1,7),WIN(2,7),WIN(3,7),WIN(4,7))
          VALUE(7)=TEXTBCK
          IF(NSTAR.GT.0) GOTO 600
          GOTO 800
       ELSE IF(CH.EQ.'X') THEN    ! Boton derecho del raton
          CALL PGSLW(1)
          CALL PGSCI(1)
          CALL PGPTXT(WIN(1,7),WIN(3,7)+0.01,0.0,0.0,VALUE(7))
          CALL PGPTXT(WIN(1,8),WIN(3,8)+0.01,0.0,0.0,VALUE(8))
          CALL PGPTXT(WIN(1,9),WIN(3,9)+0.01,0.0,0.0,VALUE(9))
          CALL PGPTXT(WIN(1,10),WIN(3,10)+0.01,0.0,0.0,VALUE(10))          
          NNOBS=NNOBS+1
          IF(NNOBS.GT.10) NNOBS=1
          VALUE(7)=OBSNAME(NNOBS)
          VALUE(8)=OBSLONG(NNOBS)
          CALL CH2DEC(VALUE(8),LONG,1)
          LONG=LONG/15.D0
          VALUE(9)=OBSLAT(NNOBS)
          CALL CH2DEC(VALUE(9),LAT,2)
          VALUE(10)=OBSHEI(NNOBS)
          READ(VALUE(10)(1:4),*)ALTOBS
          CALL PGSCI(4)
          CALL PGPTXT(WIN(1,7),WIN(3,7)+0.01,0.0,0.0,VALUE(7))
          CALL PGPTXT(WIN(1,8),WIN(3,8)+0.01,0.0,0.0,VALUE(8))
          CALL PGPTXT(WIN(1,9),WIN(3,9)+0.01,0.0,0.0,VALUE(9))
          CALL PGPTXT(WIN(1,10),WIN(3,10)+0.01,0.0,0.0,VALUE(10))          
          CALL PGSCI(1)
          LUNPLOT=.TRUE.
          GOTO 600
       ELSE
          GOTO 70
       END IF
C
C----------------------------------------------------------------------
C     Rellenamos el valor de la longitud del Observatorio
C
 800   CALL PGSCI(4)
       CALL PGSLW(8)
       CALL PGRECT(WIN(1,8),WIN(2,8),WIN(3,8),WIN(4,8))
       CALL PGRSTR(WIN(1,8),WIN(3,8)+0.01,VALUE(8),TEXTBCK,NN,0)
       CALL PGSLW(8)
       CALL PGSCI(1)
       CALL PGRECT(WIN(1,8),WIN(2,8),WIN(3,8),WIN(4,8))
       VALUE(8)=TEXTBCK
C
       IF(VALUE(8)(1:1).NE.'-'.AND.VALUE(8)(1:1).NE.'+'
     >                        .AND.VALUE(8)(1:1).NE.' ') GOTO 850
       IF(VALUE(8)(5:5).NE.' '.OR.VALUE(8)(8:8).NE.' ') GOTO 850
       DO J=2,10
          IF(TESTN(VALUE(8)(J:J)).EQ.0.AND.
     >       VALUE(8)(J:J).NE.' ') GOTO 850 
       END DO
C
       CALL CH2DEC(VALUE(8),LONG,1)
       LONG=LONG/15.D0
C
       IF(NN.EQ.7) THEN
          GOTO 900
       ELSE
          LUNPLOT=.TRUE.
          GOTO 600
       ENDIF
C
 850   CALL IOWIN1(' Please use format :  +000 00 00        ')
       GOTO 800
C----------------------------------------------------------------------
C     Rellenamos el valor de la latitud del Observatorio
C
 900   CALL PGSCI(4)
       CALL PGRECT(WIN(1,9),WIN(2,9),WIN(3,9),WIN(4,9))
       CALL PGRSTR(WIN(1,9),WIN(3,9)+0.01,VALUE(9),TEXTBCK,NN,0)
       CALL PGSCI(1)
       CALL PGSLW(8)
       CALL PGRECT(WIN(1,9),WIN(2,9),WIN(3,9),WIN(4,9))
       VALUE(9)=TEXTBCK
C
       IF(VALUE(9)(1:1).NE.'-'.AND.VALUE(9)(1:1).NE.'+'
     >                        .AND.VALUE(9)(1:1).NE.' ') GOTO 950
       IF(VALUE(9)(4:4).NE.' '.OR.VALUE(9)(7:7).NE.' ') GOTO 950
       DO J=2,9
          IF(TESTN(VALUE(9)(J:J)).EQ.0.AND.
     >             VALUE(9)(J:J).NE.' ') GOTO 950 
       END DO
C
       CALL CH2DEC(VALUE(9),LAT,2)
C
       IF(NN.EQ.7) THEN 
          GOTO 1000
       ELSE
          LUNPLOT=.TRUE.
          GOTO 600
       END IF
C
 950   CALL IOWIN1(' Please use format :  +00 00 00         ')
       GOTO 900
C----------------------------------------------------------------------
C     Rellenamos el valor de la altura del Observatorio
C
 1000  CALL PGSCI(4)
       CALL PGRECT(WIN(1,10),WIN(2,10),WIN(3,10),WIN(4,10))
       CALL PGRSTR(WIN(1,10),WIN(3,10)+0.01,VALUE(10),TEXTBCK,NN,0)
       CALL PGSCI(1)
       CALL PGSLW(8)
       CALL PGRECT(WIN(1,10),WIN(2,10),WIN(3,10),WIN(4,10))
       CALL PGSLW(1)
       VALUE(10)=TEXTBCK
C
       DO J=1,4
          IF(TESTN(VALUE(10)(J:J)).EQ.0.AND.
     >       VALUE(10)(J:J).NE.' ') GOTO 1050 
       END DO
C
       IF(VALUE(10)(1:4).EQ.'    ') THEN
          ALTOBS=0.
       ELSE
          READ(VALUE(10)(1:4),*) ALTOBS
       END IF
C
       LUNPLOT=.TRUE.
       NN=10
       GOTO 600
C
 1050  CALL IOWIN1(' Please use format :  0000m             ')
       GOTO 1000
C----------------------------------------------------------------------
C     Cambiamos la fecha
C
 1100  CALL PGSCI(4) 
       CALL PGSLW(8)
       CALL PGRECT(WIN(1,11),WIN(2,11),WIN(3,11),WIN(4,11))
       CALL PGRSTR(WIN(1,11),WIN(3,11)+0.01,VALUE(11),TEXTBCK,NN,0)
       CALL PGSCI(1)
       CALL PGSLW(8)
       CALL PGRECT(WIN(1,11),WIN(2,11),WIN(3,11),WIN(4,11))
       CALL PGSLW(1)
       VALUE(11)=TEXTBCK
C
       DO J=1,2
          IF(TESTN(VALUE(11)(J:J)).EQ.0.) GOTO 1150
       END DO
       DO J=4,5
          IF(TESTN(VALUE(11)(J:J)).EQ.0.) GOTO 1150
       END DO
       DO J=7,10
          IF(TESTN(VALUE(11)(J:J)).EQ.0.) GOTO 1150
       END DO
C
       READ(VALUE(11)(1:2),'(I2)')IDAY
       READ(VALUE(11)(4:5),'(I2)')IMONTH
       READ(VALUE(11)(7:10),'(I4)')IYEAR
C
       LUNPLOT=.TRUE.
       GOTO 600
C
 1150  CALL IOWIN1(' Please use format :  DD-MM-YYYY        ')
       GOTO 1100
C----------------------------------------------------------------------
C     Escribimos el fichero con el listado de las estrellas.
C     No Orden, Nombre, A.R., Dec. y Equinocio.
C
 1200  CALL PGSLW(12)
       CALL PGSCI(4)
       CALL PGRECT(WIN(1,12),WIN(2,12),WIN(3,12),WIN(4,12))
       OPEN(20,FILE='ves.lis',STATUS='UNKNOWN',
     >         FORM='FORMATTED')
       DO I=1,NSTAR
          WRITE(20,'(I2,A1,A11,A1,A9,A1,A9,A1,A6)')
     >     I,',',STARNAME(I),',',STARRA(I),',',STARDEC(I),',',EQUINOX(I)
       END DO
C
       CLOSE(20)
       CALL PGSCI(1)
       CALL PGRECT(WIN(1,12),WIN(2,12),WIN(3,12),WIN(4,12))
       GOTO 70
C----------------------------------------------------------------------
C     Generamos un fichero postscript con la informacion que tenemos
C     actualmente en pantalla
 1300  CALL PGSLW(12)
       CALL PGSCI(4)
       CALL PGRECT(WIN(1,13),WIN(2,13),WIN(3,13),WIN(4,13))
       CALL PGSLW(1)
       CALL PGSCI(1)
C
       DEV(2)=PGOPEN('ves.ps/vps')
       CALL PGSLCT(DEV(2))
       CALL PGSVP(0.,1.,0.,1.)
       CALL PGSWIN(0.,1.,0.,1.)
       CALL PGSCF(2)
       CALL PGSCH(0.55)
       CALL PGTEXT(0.85,0.016,COPYR)
       CALL PGSFS(2)                    ! Dibujamos un doble borde exterior
       CALL PGRECT(.01,.99,.01,.99)
       CALL PGRECT(.014,.986,.013,.987)
C Escribimos la informacion del observatorio
       CALL PGSCH(0.8)   
       CALL PGTEXT(.125,.960,'Observatory :')
       CALL PGTEXT(.125,.942,'Longitud :')
       CALL PGTEXT(.125,.924,'Latitud :')
       CALL PGTEXT(.125,.906,'Height :')
       CALL PGPTXT(.400,.960,0.0,1.0,VALUE(7))
       CALL PGPTXT(.400,.942,0.0,1.0,VALUE(8))
       CALL PGPTXT(.400,.924,0.0,1.0,VALUE(9))
       CALL PGPTXT(.400,.906,0.0,1.0,VALUE(10))
C Escribimos info de la fecha
       CALL PGTEXT(.600,.960,'Date :')  
       CALL PGPTXT(.875,.960,0.0,1.0,VALUE(11))
       CALL PGTEXT(.600,.942,'JD (0h UT) :')
       WRITE(CFJ0,'(F10.2)')FJ0
       CALL PGPTXT(.875,.942,0.0,1.0,CFJ0)
       CALL PGTEXT(.600,.924,'Night length :')
       CALL PGTEXT(.600,.906,'Twilight free :')
       NIGHTL=24.-TUOCASO2+TUORTO2
       IF(NIGHTL.GT.24.) NIGHTL=NIGHTL-24.
       NIGHTLF=24.-TUCREPV2+TUCREPM2
       IF(NIGHTLF.GT.24.) NIGHTLF=NIGHTLF-24.
       CALL R2CH(NIGHTL,CNIGHTL)
       CALL R2CH(NIGHTLF,CNIGHTLF)
       CALL PGPTXT(.875,.924,0.0,1.0,CNIGHTL)
       CALL PGPTXT(.875,.906,0.0,1.0,CNIGHTLF)
       CALL PGMOVE(.025,.900)
       CALL PGDRAW(.975,.900)
       CALL PGMOVE(.025,.897)
       CALL PGDRAW(.975,.897)
C Dibujamos el plot.
       CALL PGSVP(.075,.925,.340,.840)
       CALL PGSWIN(XMIN,XMAX,YMIN,YMAX)
       CALL PGSCH(0.6)
       CALL PGPTXT(XMIN,91.,0.0,0.0,'Sunset')
       CALL R2CH(TUOCASO2,CHOCASOTU)
       CALL PGPTXT(XMIN,93.,0.0,0.0,CHOCASOTU//' UT')
       CALL R2CH(TSLOCASO,CHOCASOTSL)
       CALL PGPTXT(XMIN,95.,0.0,0.0,CHOCASOTSL//' LST')
       CALL PGPTXT(XMAX,91.,0.0,1.0,'Sunrise')
       CALL R2CH(TUORTO2,CHORTOTU)
       CALL PGPTXT(XMAX,93.,0.0,1.0,CHORTOTU//'  UT ')
       CALL R2CH(TSLORTO,CHORTOTSL)
       CALL PGPTXT(XMAX,95.,0.0,1.0,CHORTOTSL//' LST')
       CALL PGSCH(1.0)
C
       IF(LUT) THEN  ! Altura frente a T.U.
          CALL PGTBOX('ABCNTSZYXHO',0.0,0,'ABCN',0.0,0)
          CALL PGPTXT(XMIN-2200.,42.,90.,0.0,'Elevation')
          CALL PGPTXT(XMIN-2000.,-4.5,0.,0.0,'U.T.')
          CALL PGSLS(4) 
          DO I=10,90,10      ! lineas horizontales cada 10 grados de elevac.
             CALL PGMOVE(XMIN,REAL(I))
             CALL PGDRAW(XMAX,REAL(I))
          END DO
          DO I=0,172800,3600        ! Lineas de puntos verticales cada hora
             IF(I.GT.XMIN.AND.I.LT.XMAX) THEN
                CALL PGMOVE(REAL(I),0.)
                CALL PGDRAW(REAL(I),90.)
             END IF
          END DO
          CALL PGSLS(2)  ! Lineas verticales a trazos para marcar los crepusc
          CALL PGMOVE(XCREPV,0.)
          CALL PGDRAW(XCREPV,90.)
          CALL PGMOVE(XCREPM,0.)
          CALL PGDRAW(XCREPM,90.)
          CALL PGSLS(1)
          DO I=1,10   ! Ticks para marcar las sec Z (masas de aire)  
             Y=ACOS(1./SECZ(I))
             Y=Y*180./PI         
             Y=90.-Y
             CALL PGMOVE(XMAX-900.,Y)
             CALL PGDRAW(XMAX,Y)
             CALL PGPTXT(XMAX+300.,Y-1.,0.0,0.0,CHSECZ(I))
          END DO
          CALL PGPTXT(XMAX+2000.,72.,90.,0.0,'sec(Z)')
          DO J=1,NSTAR  ! Pintamos las trazas de las estrellas
             DO I=1,100
                TUPSU(I)=TUPS(J,I)
                ALTU(I)=ALT(J,I)
             END DO
             CALL PGSLW(3)
             CALL PGLINE(100,TUPSU,ALTU)
             CALL PGSLW(1)
             WRITE(CH2,'(I2)') J
             CALL PGSCH(.5)
             DO I=1,9
                IF(ALTU(10*I+J).GT.0.AND.ALTU(10*I+J).LT.89) THEN
                   CALL PGTEXT(TUPSU(10*I+J),ALTU(10*I+J),CH2)
                END IF
             END DO
             CALL PGSCH(1.)
          END DO
C     Pintamos la Luna si se ve
          DO I=1,48
             IF(LMOON(I)) THEN 
                CALL PGSCH(1.5)
                YYY=MOONALT(I)
                CALL PGPT1(I*3600.,YYY,2291)
                WRITE(CHF,'(F4.2)')FASE(I)
                CALL PGSCH(0.6)
                CALL PGTEXT(I*3600.,YYY-0.5,CHF)
             END IF
          END DO
       ELSE IF(LAZI) THEN
          CALL PGSVP(.100,.900,.325,.875)
          CALL PGSWIN(-1.,1.,-1.,1.)
          CALL PGCIRC(0.,0.,0.9)
          CALL PGSLS(2)  ! Trazamos cruz con lineas de guiones
          CALL PGMOVE(-.9,0.) 
          CALL PGDRAW(.9,0.)
          CALL PGMOVE(0.,-0.9)
          CALL PGDRAW(0.,0.9)
          CALL PGSCH(2.)
          CALL PGTEXT(-.05,0.92,'S')
          CALL PGTEXT(-.05,-1.01,'N')
          CALL PGTEXT(-1.,0.,'E')
          CALL PGTEXT(0.92,0.,'W')
          CALL PGSLS(4)
          DO I=1,5 ! Circulos de puntos cada 15 grados
             CALL PGCIRC(0.,0.,REAL(I)*0.15)
          END DO 
          CALL PGSCH(0.7)
          CALL PGTEXT(0.0,-.74,'15')
          CALL PGTEXT(0.0,-.59,'30')
          CALL PGTEXT(0.0,-.44,'45')
          CALL PGTEXT(0.0,-.29,'60')
          CALL PGTEXT(0.0,-.14,'75')
          CALL PGSCH(1.)
          CALL PGSLS(1)
          DO J=1,NSTAR
             WRITE(CH2,'(I2)') J
             CALL PGSLW(3)
             DO I=1,100
                RADIO(I)=SQRT(XAZI(J,I)**2+YAZI(J,I)**2)
                IF(I.GT.1.AND.RADIO(I).LT.0.9.AND.RADIO(I).GT.0.)THEN
                   CALL PGMOVE(XAZI(J,I-1),YAZI(J,I-1))
                   CALL PGDRAW(XAZI(J,I),YAZI(J,I))
                END IF
             END DO
             CALL PGSLW(1)
             CALL PGSCH(0.5)
             DO I=1,9
                IF(RADIO(10*I).LT.0.9.AND.RADIO(10*I).GT.0.) THEN
                   CALL PGTEXT(XAZI(J,10*I),YAZI(J,10*I),CH2)
                END IF
             END DO
             CALL PGSCH(1.0)
          END DO
C     Pintamos la Luna si se ve
          DO I=1,48
             IF(LMOON(I)) THEN
                CALL PGSCH(1.5)
                CALL PGPT1(XMAZI(I),YMAZI(I),2291)
                WRITE(CHF,'(F4.2)')FASE(I)
                CALL PGSCH(0.6)
                CALL PGTEXT(XMAZI(I),YMAZI(I)-0.015,CHF)
             END IF
          END DO
       END IF
C
       CALL PGSVP(0.,1.,0.,1.)
       CALL PGSWIN(0.,1.,0.,1.)
       CALL PGSLS(1)
       CALL PGMOVE(.025,.300)
       CALL PGDRAW(.975,.300)
       CALL PGMOVE(.025,.297)
       CALL PGDRAW(.975,.297)
C     Escribimos la lista de objetos 
       IF(NSTAR.GT.0) THEN
          CALL PGSCH(.6)
          CALL PGTEXT(.050,.280,'  Object')
          CALL PGTEXT(.235,.280,'\ga   Initial    \gd  ')
          CALL PGTEXT(.450,.280,'\ga   '//CTFF//'   \gd  ')
          CALL PGTEXT(.640,.280,'Culmination  U.T.')
          CALL PGTEXT(.840,.280,'Visible')
          CALL PGMOVE(.050,.275)
          CALL PGDRAW(.950,.275)
          DO I=1,NSTAR
             IF(I.EQ.17) THEN
                CALL IOWIN1(' More than 16 objects                    ')
                GOTO 1350
             END IF
             WRITE(CH2,'(I2)') I
             CALL PGPTXT(.050,.275-I*.015,0.,0.,CH2)
             CALL PGPTXT(.070,.275-I*.015,0.,0.,STARNAME(I))
             CALL PGPTXT(.280,.275-I*.015,0.,1.,STARRA(I))
             CALL PGPTXT(.375,.275-I*.015,0.,1.,STARDEC(I))
             CALL RR2CH(ARFF(I),STARRAF(I))
             CALL PGPTXT(.500,.275-I*.015,0.,1.,STARRAF(I))
             CALL RR2CH(DECFF(I),STARDECF(I))
             CALL PGPTXT(.600,.275-I*.015,0.,1.,STARDECF(I))
C   Calculamos la altura, direccion y T.U. de la culminacion.
             CALL CULMINA(LONG,LAT,TS0,ARFF(I),DECFF(I),
     >                    ALTCUL,DIRCUL,TUCUL)
             IF(DIRCUL.NE.'I') THEN
                WRITE(CALTCUL,'(F5.1)') ALTCUL
                CALL PGPTXT(.640,.275-I*.015,0.,0.,CALTCUL)
                CALL PGPTXT(.685,.275-I*.015,0.,0.,DIRCUL)
                CALL R2CH(TUCUL,CTUCUL)
                CALL PGPTXT(.785,.275-I*.015,0.,1.,CTUCUL)
             END IF
C   Calculamos la visibilidad del objeto
             CALL ORCASO(DECFF(I),90.D0,ESTADO)
             IF(ESTADO.EQ.'I') THEN
                CALL PGTEXT(.835,.275-I*.015,'Invisible')
             ELSE IF(ESTADO.EQ.'C') THEN
                CALL PGTEXT(.835,.275-I*.015,'Circumpolar')
             ELSE
                CALL TIEMPOU(ARFF(I),HORTO,LONG,TS0,TSL,TUORSTAR)
                CALL TIEMPOU(ARFF(I),HOCASO,LONG,TS0,TSL,TUOCSTAR)
                IF(TUOCSTAR.GT.TUORSTAR) THEN
                   IF((TUORSTAR+24.)*3600.LT.TVISFIN(I)) THEN
                      TUORSTAR=TUORSTAR+24.
                      TUOCSTAR=TUOCSTAR+24.
                   END IF
                   TVISINI(I)=DMAX1(TVISINI(I),TUORSTAR*3600.D0)
                   TVISFIN(I)=DMIN1(TVISFIN(I),TUOCSTAR*3600.D0)
                ELSE
                   TUOCSTAR=TUOCSTAR+24.
                   TVISINI(I)=DMAX1(TVISINI(I),TUORSTAR*3600.D0)
                   TVISFIN(I)=DMIN1(TVISFIN(I),TUOCSTAR*3600.D0)
                END IF
                TVIS(I)=(TVISFIN(I)-TVISINI(I))/3600.D0
                CALL R2CH(TVIS(I),CTVIS)
                CALL PGPTXT(.900,.275-I*.015,0.,1.,CTVIS)
             END IF
          END DO
       END IF
C
 1350  CALL PGSLS(1)
       CALL PGSCH(1.)
C
       CALL PGCLOS(DEV(2))
       CALL PGSLCT(DEV(1))
       CALL PGSLW(12)
       CALL PGRECT(WIN(1,13),WIN(2,13),WIN(3,13),WIN(4,13))
       CALL PGSLW(1)
       GOTO 70

C----------------------------------------------------------------------
C     Salimos del programa.
C
 1400  CALL PGSLW(12)
       CALL PGSCI(4)
       CALL PGSLS(1)
       CALL PGRECT(WIN(1,14),WIN(2,14),WIN(3,14),WIN(4,14))
       CALL PGSLW(1)
       CALL PGSCI(1)
       GOTO 2000
C
C----------------------------------------------------------------------
C     Selecionamos la altura sobre el horizonte como coordenada Y.
C
 1500  CALL PGSCI(1)
       CALL PGSLW(12)
       CALL PGSLS(1)
       CALL PGRECT(WIN(1,16),WIN(2,16),WIN(3,16),WIN(4,16))
       CALL PGSCI(4)
       CALL PGRECT(WIN(1,15),WIN(2,15),WIN(3,15),WIN(4,15))
       CALL PGSLW(1)
C
       LSECZ=.FALSE.
       LELEV=.TRUE.
C
       IF(LUT) THEN  
          GOTO 60
       ELSE   ! Azimut
          CALL PGSLS(4)
          CALL PGSVP(0.24,0.82,0.20,0.94)
          CALL PGSWIN(-1.,1.,-1.,1.)
          CALL PGSCI(0)   ! Borramos los circulos de puntos que marcan las secZ
          CALL PGSCH(.7)
          DO I=1,10              
             Y=ACOS(1./SECZ(I))
             Y=Y*1.8/PI         
             IF(I.EQ.1) CALL PGTEXT(0.,Y,'1.0')
             IF(I.EQ.2) CALL PGTEXT(0.,Y,'1.1')
             IF(I.EQ.7) CALL PGTEXT(0.,Y,'2.0')
             CALL PGCIRC(0.,0.,Y)
          END DO
          CALL PGSCI(1)  ! Pintamos los circulos de puntos de altura
          CALL PGSCH(1.)
          DO I=1,5
             CALL PGCIRC(0.,0.,REAL(I)*0.150)
          END DO 
          CALL PGSCH(0.7)
          CALL PGTEXT(0.0,-.74,'15')
          CALL PGTEXT(0.0,-.59,'30')
          CALL PGTEXT(0.0,-.44,'45')
          CALL PGTEXT(0.0,-.29,'60')
          CALL PGTEXT(0.0,-.14,'75')
          CALL PGSCH(1.)
          CALL PGSVP(0.0,1.0,0.0,1.0)
          CALL PGSWIN(0.0,2.0,0.0,1.0)
       END IF
C
       GOTO 70
C----------------------------------------------------------------------
C     Selecionamos sec Z (aprox. masa de aire) como coordenada Y. 
C
 1600  CALL PGSCI(1)
       CALL PGSLS(1)
       CALL PGSLW(12)
       CALL PGRECT(WIN(1,15),WIN(2,15),WIN(3,15),WIN(4,15))
       CALL PGSCI(4)
       CALL PGRECT(WIN(1,16),WIN(2,16),WIN(3,16),WIN(4,16))
       CALL PGSLW(1)
C
       LELEV=.FALSE.
       LSECZ=.TRUE. 
C
       IF(LUT) THEN  
          GOTO 60
       ELSE   ! Azimut
          CALL PGSLS(1)
          CALL PGSVP(0.24,0.82,0.20,0.94)
          CALL PGSWIN(-1.,1.,-1.,1.)
          CALL PGSCI(0)  ! Borramos los circulos de puntos de altura
          CALL PGSLW(10)
          CALL PGSCH(1.)
          DO I=1,5
             CALL PGCIRC(0.,0.,REAL(I)*0.150)
          END DO 
          CALL PGSLW(1)
          CALL PGSCH(0.7)
          CALL PGTEXT(0.0,-.74,'15')
          CALL PGTEXT(0.0,-.59,'30')
          CALL PGTEXT(0.0,-.44,'45')
          CALL PGTEXT(0.0,-.29,'60')
          CALL PGTEXT(0.0,-.14,'75')
          CALL PGSLS(4)
          CALL PGSCI(1)   ! Pintamos los circulos de puntos que marcan las secZ
          DO I=1,10              
             Y=ACOS(1./SECZ(I)) ! Circulos  para marcar las
             Y=Y*1.8/PI         ! sec Z (masas de aire)
             IF(I.EQ.1) CALL PGTEXT(0.,Y,'1.0')
             IF(I.EQ.2) CALL PGTEXT(0.,Y,'1.1')
             IF(I.EQ.7) CALL PGTEXT(0.,Y,'2.0')
             CALL PGCIRC(0.,0.,Y)
          END DO
          CALL PGSCH(1.)
          CALL PGSVP(0.0,1.0,0.0,1.0)
          CALL PGSWIN(0.0,2.0,0.0,1.0)
       END IF
C
       GOTO 70
C----------------------------------------------------------------------
C     Selecionamos Tiempo Universal como coordenada X
C
 1700  CALL PGSCI(1)
       CALL PGSLW(12)
       CALL PGSLS(1)
       CALL PGRECT(WIN(1,18),WIN(2,18),WIN(3,18),WIN(4,18))
       CALL PGSCI(4)
       CALL PGRECT(WIN(1,17),WIN(2,17),WIN(3,17),WIN(4,17))
       CALL PGSLW(1)
C
       LAZI=.FALSE.
       LUT=.TRUE.
C
       CALL PGSCI(0)  ! Borramos el grafico de azimut
       CALL PGSVP(0.24,0.82,0.20,0.94)
       CALL PGSWIN(-1.,1.,-1.,1.)
       CALL PGSLS(1)
       CALL PGCIRC(0.,0.,0.9)
       CALL PGSLS(2)
       CALL PGMOVE(-.9,0.)
       CALL PGDRAW(.9,0.)
       CALL PGMOVE(0.,-0.9)
       CALL PGDRAW(0.,0.9)
       CALL PGSCH(2.)
       CALL PGTEXT(-.05,0.92,'S')
       CALL PGTEXT(-.05,-1.01,'N')
       CALL PGTEXT(-1.,0.,'E')
       CALL PGTEXT(0.92,0.,'W')
       CALL PGSCH(1.)
C
       CALL PGSLS(4)
       CALL PGSLW(1)
       IF(LELEV) THEN
          DO I=1,5
             CALL PGCIRC(0.,0.,REAL(I)*0.15)
          END DO 
          CALL PGSCH(0.7)
          CALL PGTEXT(0.0,-.74,'15')
          CALL PGTEXT(0.0,-.59,'30')
          CALL PGTEXT(0.0,-.44,'45')
          CALL PGTEXT(0.0,-.29,'60')
          CALL PGTEXT(0.0,-.14,'75')
          CALL PGSCH(1.)
       ELSE 
          CALL PGSCH(0.7)
          DO I=1,10              
             Y=ACOS(1./SECZ(I)) ! Lineas horizontales para marcar las
             Y=Y*1.8/PI         ! sec Z (masas de aire)
             IF(I.EQ.1) CALL PGTEXT(0.,Y,'1.0')
             IF(I.EQ.2) CALL PGTEXT(0.,Y,'1.1')
             IF(I.EQ.7) CALL PGTEXT(0.,Y,'2.0')
             CALL PGCIRC(0.,0.,Y)
          END DO
          CALL PGSCH(1.)
       END IF
C     Borramos la Luna
       DO I=1,48
          IF(LMOON(I)) THEN
             CALL PGSCH(1.5)
             CALL PGPT1(XMAZI(I),YMAZI(I),2291)
             WRITE(CHF,'(F4.2)')FASE(I)
             CALL PGSCH(0.6)
             CALL PGTEXT(XMAZI(I),YMAZI(I)-0.01,CHF)
          END IF
       END DO
       CALL PGSCH(1.)
C
       NSTAROLD=NSTAR
       NSTAR=0
C
       CALL PGSCI(1)
       CALL PGSLW(1)
c
       GOTO 60
C----------------------------------------------------------------------
C     Selecionamos Azimut como coordenada X
C
 1800  CALL PGSCI(1)
       CALL PGSLW(12)
       CALL PGRECT(WIN(1,17),WIN(2,17),WIN(3,17),WIN(4,17))
       CALL PGSCI(4)
       CALL PGRECT(WIN(1,18),WIN(2,18),WIN(3,18),WIN(4,18))
       CALL PGSLW(1)
C
       LUT=.FALSE.
       LAZI=.TRUE.
C
       CALL PGSCI(0)         ! Borramos
       CALL PGSVP(0.22,0.85,0.22,0.90)
       CALL PGSWIN(XMIN,XMAX,YMIN,YMAX)
C
       IF(LSECZ) THEN
          CALL PGSLS(1) 
          CALL PGTBOX('ABCNTSZYXHO',0.0,0,'ABC',0.0,0)
          CALL PGSLS(4) 
          DO I=1,10              
             Y=ACOS(1./SECZ(I)) ! Lineas horizontales para marcar las
             Y=Y*180./PI        ! sec Z (masas de aire)
             Y=90.-Y
             CALL PGMOVE(XMIN,Y)
             CALL PGDRAW(XMAX,Y)
             CALL PGPTXT(XMIN-300.,Y-1.,0.0,1.0,CHSECZ(I))
          END DO
       ELSE IF(LELEV) THEN
          CALL PGSLS(1) 
          CALL PGTBOX('ABCNTSZYXHO',0.0,0,'ABCNTS',0.0,0)
          CALL PGSLS(4) 
            DO I=10,90,10      ! lineas horizontales cada 10 grados de elevac.
             CALL PGMOVE(XMIN,REAL(I))
             CALL PGDRAW(XMAX,REAL(I))
          END DO
       END IF
C
       CALL PGSLS(4)
       DO I=0,172800,3600        ! Lineas de puntos verticales cada hora
          IF(I.GT.XMIN.AND.I.LT.XMAX) THEN
             CALL PGMOVE(REAL(I),0.)
             CALL PGDRAW(REAL(I),90.)
          END IF
       END DO
C
       CALL PGSLS(2)  ! Lineas verticales a trazos para marcar los crepusc
       CALL PGMOVE(XCREPV,0.)
       CALL PGDRAW(XCREPV,90.)
       CALL PGMOVE(XCREPM,0.)
       CALL PGDRAW(XCREPM,90.)
C
       CALL PGSCI(1) ! Pintamos la trama
       CALL PGSVP(0.24,0.82,0.20,0.94)
       CALL PGSWIN(-1.,1.,-1.,1.)
       CALL PGSLS(1)
       CALL PGCIRC(0.,0.,0.9)
       CALL PGSLS(2)  ! Trazamos cruz con lineas de guiones
       CALL PGMOVE(-.9,0.)
       CALL PGDRAW(.9,0.)
       CALL PGMOVE(0.,-0.9)
       CALL PGDRAW(0.,0.9)
       CALL PGSCH(2.)
       CALL PGTEXT(-.05,0.92,'S')
       CALL PGTEXT(-.05,-1.01,'N')
       CALL PGTEXT(-1.,0.,'E')
       CALL PGTEXT(0.92,0.,'W')
C
       CALL PGSLS(4) 
       CALL PGSLW(1)
       CALL PGSCH(1.)
       IF(LELEV)THEN
          DO I=1,5 ! Circulos de puntos cada 15 grados
             CALL PGCIRC(0.,0.,REAL(I)*0.15)
          END DO 
          CALL PGSCH(0.7)
          CALL PGTEXT(0.0,-.74,'15')
          CALL PGTEXT(0.0,-.59,'30')
          CALL PGTEXT(0.0,-.44,'45')
          CALL PGTEXT(0.0,-.29,'60')
          CALL PGTEXT(0.0,-.14,'75')
          CALL PGSCH(1.)
       ELSE 
          CALL PGSCH(0.7)
          DO I=1,10              
             Y=ACOS(1./SECZ(I)) ! Circulos de puntos  para marcar las
             Y=Y*1.8/PI         ! sec Z (masas de aire)
             IF(I.EQ.1) CALL PGTEXT(0.,Y,'1.0')
             IF(I.EQ.2) CALL PGTEXT(0.,Y,'1.1')
             IF(I.EQ.7) CALL PGTEXT(0.,Y,'2.0')
             CALL PGCIRC(0.,0.,Y)
          END DO
          CALL PGSCH(1.0)
       END IF
C
       CALL PGSCI(4)  !Pintamos las trazas de los objetos
       CALL PGSLS(1)
       CALL PGSLW(2)
       DO J=1,NSTAR
          WRITE(CH2,'(I2)') J
          CALL PGSVP(0.24,0.82,0.20,0.94)
          CALL PGSWIN(-1.,1.,-1.,1.)
          DO I=1,100
             RADIO(I)=SQRT(XAZI(J,I)**2+YAZI(J,I)**2)
             IF(I.GT.1.AND.RADIO(I).LT.0.9.AND.RADIO(I).GT.0.)THEN
                IF(ABS(XAZI(J,I-1)*YAZI(J,I-1)).GT.0) THEN
                   CALL PGMOVE(XAZI(J,I-1),YAZI(J,I-1))
                ELSE
                   CALL PGMOVE(XAZI(J,I),YAZI(J,I))
                END IF
                CALL PGDRAW(XAZI(J,I),YAZI(J,I))
             END IF
          END DO
          CALL PGSLW(1)
          CALL PGSCH(0.5)
          DO I=1,9
             IF(RADIO(10*I+J).LT.0.9.AND.RADIO(10*I+J).GT.0.) THEN
                CALL PGTEXT(XAZI(J,10*I+J),YAZI(J,10*I+J),CH2)
             END IF
          END DO
          CALL PGSCH(1.0)
       END DO
C      Ponemos la Luna
       CALL PGSCI(5)
       DO I=1,48
          IF(LMOON(I)) THEN
             CALL PGSCH(1.5)
             CALL PGPT1(XMAZI(I),YMAZI(I),2291)
             WRITE(CHF,'(F4.2)')FASE(I)
             CALL PGSCH(0.6)
             CALL PGTEXT(XMAZI(I),YMAZI(I)-0.01,CHF)
          END IF
       END DO
C
       CALL PGSCH(1.)
       CALL PGSCI(1)
       CALL PGSLW(1)
       CALL PGSVP(0.0,1.0,0.0,1.0)
       CALL PGSWIN(0.0,2.0,0.0,1.0)
C
       GOTO 70
C----------------------------------------------------------------------
C
 2000  CALL PGASK(.FALSE.)
       CALL PGEND
C
       STOP
       END
C
C**********************************************************************
C     Para modificar el texto que haya en una ventana.
C
      SUBROUTINE PGRSTR(X, Y, TEXT,TEXTBCK,NN,SP)
C
      IMPLICIT NONE
      INTEGER NN            ! Numero de ventana
      INTEGER LSTR          ! Longitud del texto.
      INTEGER SP            ! SP=1 pone un espacio
      REAL X, Y             ! Posicion donde empieza el texto.
      REAL XBOX(4),YBOX(4)  ! Coordenadas de la caja que enmarca el texto. 
      REAL XX,YY
      CHARACTER*1   CH
      CHARACTER*(*) TEXT
      CHARACTER*11  TEXTOLD
      CHARACTER*11  TEXTBCK
C
      TEXTOLD=TEXT
      LSTR=0
      CALL PGSLW(1)
      CALL PGSCI(1)            ! Borramos el texto viejo
      CALL PGSCH(1.)
      CALL PGPTXT (X,Y,0.0,0.0,TEXT)
      CALL PGSCI(4)
C
      IF(SP.EQ.1) THEN
         LSTR=1
         TEXT=' '
      END IF
C
 10   IF(LSTR.GT.0) THEN
         CALL PGPTXT(X,Y,0.0,0.0,TEXT(1:LSTR))
         CALL PGQTXT(X,Y,0.0,0.0,TEXT,XBOX,YBOX)
         XX=XBOX(4)
         YY=YBOX(4)
      ELSE
         XX=X
         YY=Y
      END IF
C
      CALL PGBAND(0,0,XX,YY,XX,YY,CH)
C
      CALL PGSCI(1)            ! Borramos el texto viejo
      IF(LSTR.GT.0) CALL PGPTXT (X,Y,0.0,0.0,TEXT(1:LSTR))
      CALL PGSCI(4)
C
      IF(CH.EQ.CHAR(92)) CH = '*'  ! Para evitar problemas con "\" y PGPLOT 
C
      IF(ICHAR(CH).EQ.8.) THEN  ! Backspace o ctrl H borran ultimo caracter
         IF(LSTR.GT.0) THEN                        
            TEXT(LSTR:LSTR)= ' '
            LSTR=LSTR-1
         END IF
      ELSE IF (ICHAR(CH).EQ.127) THEN  ! Del borra todo el texto
         TEXT(1:LSTR)=' '
         LSTR=0
         IF(NN.EQ.2) LSTR=1
      ELSE IF (ICHAR(CH).LT.32) THEN
         IF(LSTR.GT.0) THEN
           IF(NN.EQ.2.AND.TEXT(3:3).EQ.' ') THEN
             CALL PGPTXT(X,Y,0.0,0.0,TEXTOLD)
             TEXTBCK=TEXTOLD
             CALL IOWIN1(' Please use format :  00 00 00          ')
           ELSE
             CALL PGPTXT(X,Y,0.0,0.0,TEXT(1:LSTR))
             TEXTBCK=TEXT(1:LSTR)
           END IF
         ELSE
           CALL PGPTXT(X,Y,0.0,0.0,TEXTOLD)
           TEXTBCK=TEXTOLD
         END IF
         GOTO 20
      ELSE IF (LSTR.LT.LEN(TEXT)) THEN  ! Aniadimos un caracter al texto
         LSTR=LSTR+1
         TEXT(LSTR:LSTR)=CH
      END IF
C      
      GOTO 10
C
 20   RETURN
      END
C
C**********************************************************************
      INTEGER FUNCTION TESTN(CH)
C     Funcion que chequea si un caracter es un numero del 0 al 9.
C
      CHARACTER*1 CH
C
      IF(ICHAR(CH).GT.47.AND.ICHAR(CH).LT.58) THEN
         TESTN = 1
      ELSE
         TESTN=0
      END IF
C
      RETURN
      END
C
C **********************************************************************
C                                                      
       SUBROUTINE FJ_TS0(IYEAR,IMONTH,IDAY,RHOUR,FJ,TS0)
C
C Calcula la fecha juliana y el Tiempo Sidereo en Greenwich a 0h UT a
C partir de los datos de una fecha concreta. En la fecha juliana se
C tiene en cuenta la fraccion del dia correspondiente a la hora.
C Los calculos se realizan en doble precision para no perder informacion
C por redondeo.
C
       IMPLICIT NONE
C
       INTEGER IYEAR,IMONTH,IDAY        ! Entrada
       REAL*8 RHOUR                       ! Entrada
       REAL*8 TS0                         ! Salida
       REAL*8 FJ                        ! Salida
C
       INTEGER*4 FECHA
       REAL*8 M,A,B,Y,DT,DTS0
       REAL*8 DYEAR,DMONTH,DDAY
C
       DYEAR=DBLE(IYEAR)
       DMONTH=DBLE(IMONTH)
       DDAY=DBLE(IDAY)
       FECHA=IYEAR*10000+IMONTH*100+IDAY
       IF(FECHA.GE.15821015)THEN
         A=DINT(DYEAR/1.D2)
         B=2.D0-A+DINT(A/4.D0)
       ELSE
         B=0.D0
       END IF
       IF(IYEAR.GE.0)THEN
         A=0.D0
       ELSE
         A=-.75D0
       END IF
       IF(IMONTH.GE.3)THEN
         Y=DYEAR
         M=DMONTH
       ELSE
         Y=DYEAR-1.D0
         M=DMONTH+1.2D1
       END IF
       FJ=DINT(365.25D0*Y+A)+DINT(30.6001D0*(M+1.D0))+DDAY+B+
     >        1720994.5D0
       DT=(FJ-2451545.D0)/36525.D0
       DTS0=6.D0*3.6D3+41.D0*6.D1+50.54851D0+8640184.812866D0*DT
     >      +.093104D0*DT*DT-6.2D-6*DT*DT*DT
       DTS0=DMOD(DTS0,8.64D4)
       DTS0=DTS0/8.64D4
       DTS0=DTS0*2.4D1
       TS0=REAL(DTS0)
       IF(TS0.LT.0.)TS0=TS0+24.
       IF(TS0.GT.24.)TS0=TS0-24.
       FJ=FJ+DBLE(RHOUR/24.)
       END
C
C **********************************************************************
C                                                     
       SUBROUTINE POSSOL(FJ,ARSUN,DECSUN)
C      Calculamos AR y DEC del Sol para la FJ dada
C	   Astronomical Formulae for Calculators
C	   J. Meeus Cap. 18, pag. 79 y siguientes
C
       IMPLICIT NONE
C
       REAL*8 FJ              ! Entrada
       REAL*8 ARSUN,DECSUN    ! Salida
       REAL*8 PI,RAG
       REAL*8 ALAMBDA,PAR
       REAL*8 AR,DEC
C
       PI=3.141592654D0
       RAG=180.0D0/PI
       CALL PSOL(FJ,ALAMBDA,PAR)
       CALL ECUAT(FJ,ALAMBDA,0.D0,AR,DEC)
       ARSUN=AR*RAG/15.D0
       DECSUN=DEC*RAG
       END
C ********************************************************************
C                                                     
       SUBROUTINE PSOL(FJ,ALON,PAR)
C
C	      AL : Long. media
C	      AM : Anomalia media
C	      EXC: Excentricidad
C	      EC : Ecuacion del centro
C	      V  : Anomalia verdadera
C	      R  : Radiovector
C
       IMPLICIT DOUBLE PRECISION (A-H,O-Z)
       PI=3.141592654D0
       RAG=180.0D0/PI
       DPI=2.0D0*PI
       T=(FJ-2415020.D0)/36525.D0
C
	  AL=279.69668D0+36000.76892D0*T+.0003025D0*T*T
	  AM=358.47583D0+35999.04975D0*T-0.000150D0*T*T
     >	     +.0000033D0*T*T*T
	  AL=DMOD(AL,360.D0)/RAG
	  AM=DMOD(AM,360.D0)/RAG
	  EXC=.016751014D0-.0000418D0*T-.000000126D0*T*T
	  C1=1.919460D0-.004789D0*T-.000014D0*T*T
	  C2=.020094D0-.000100D0*T
	  EC=(C1*DSIN(AM)+C2*DSIN(2.*AM)+.000293D0*DSIN(3.*AM))/RAG
C
	  ALON=AL+EC
	  V=AM+EC
	  R=1.0000002D0*(1.D0-EXC*EXC)/(1.D0+EXC*DCOS(V))
C
C   Perturbaciones
	  A=(153.23D0+22518.7541D0*T)/RAG
	  B=(216.57D0+45037.5082D0*T)/RAG
	  C=(312.69D0+32964.3577D0*T)/RAG
	  D=(350.74D0+445267.1142D0*T-.00144D0*T*T)/RAG
	  E=(231.19D0+20.20D0*T)/RAG
	  H=(353.40D0+65928.7155D0*T)/RAG
C
	  A=DMOD(A,DPI)
	  B=DMOD(B,DPI)
	  C=DMOD(C,DPI)
	  D=DMOD(D,DPI)
	  H=DMOD(H,DPI)
	  E=DMOD(E,DPI)
C
	  ALON=ALON+(.00134D0*DCOS(A)+.00154D0*DCOS(B)
     >		   +.00200D0*DCOS(C)+.00179D0*DSIN(D)
     >		   +.00178D0*DSIN(E))/RAG
	  R=R+.00000543D0*DSIN(A)+.00001575D0*DSIN(B)
     > 	     +.00001627D0*DSIN(C)+.00003076D0*DCOS(D)
     >	     +.00000927D0*DSIN(H)
C
c	  IF(IAP .EQ. 1) THEN                 ! Por si se necesitan
c	     OM=(259.18D0-1934.142D0*T)/RAG   ! coordenadas aparentes.
c	     OM=DMOD(OM,DPI)
c	     IF(OM .LT. 0.) OM=OM+DPI
c	     ALON=ALON-(0.00569D0+0.00479D0*DSIN(OM))/RAG
c	  ENDIF
	  RR=R*1.49597870D11/6.37814D6
	  PAR=ASIN(1.D0/RR)*RAG
C
       RETURN
       END
C
C
C
C ********************************************************************
C      Calculamos Ascension Recta y Declinacion       
C
       SUBROUTINE ECUAT(FJ,AL,B,A,D)
C
       IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
       DPI=6.283185307D0
       CALL EPSIL(FJ,EPSI)
C
       X0=DCOS(B)*DCOS(AL)
       Y0=DCOS(B)*DSIN(AL)
       Z0=DSIN(B)
C
       X=X0
       Y=Y0*DCOS(EPSI)-Z0*DSIN(EPSI)
       Z=Y0*DSIN(EPSI)+Z0*DCOS(EPSI)
C
       A=ATAN2(Y,X)
       IF(A .LT. 0.) A=A+DPI
       D=ASIN(Z)
C
       RETURN
       END
C**********************************************************************      
C      Calculo de la oblicuidad de la ecliptica
C
       SUBROUTINE EPSIL(FJ,EPSI)
C
       IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
       U=(FJ-2451545.D0)/3652500.D0
C
       A1=2.18D0-3375.70D0*U+0.36D0*U**2
       A2=3.51D0+125666.39D0*U+0.10D0*U**2
       EPSI=0.4090928D0-0.0226938D0*U-75.D-7*U**2+96926.D-7*U**3
     >	 -2491.D-7*U**4-12104.D-7*U**5+(446*DCOS(A1)+28*DCOS(A2))*1.D-7
       RETURN
       END
C
C **********************************************************************
C     Calculamos para un objeto (sol) el orto y ocaso, condicion en la
C     cual obtiene una distancia cenital igual a DISZEN.Si el objeto
C     es circumpolar estado=C, si es invisible estado=I y sino estado=O
C                                                      
       SUBROUTINE ORCASO(DECSUN,DISZEN,ESTADO)
C
       IMPLICIT NONE
C
       REAL*8 DECSUN,DISZEN   ! Entrada
       REAL*8 LAT,LONG        ! Entrada
       REAL*8 HORTO,HOCASO    ! Salida
       REAL*8 AORTO,AOCASO    ! Salida
C
       REAL*8 PI
       PARAMETER(PI=3.141592654D0)
       REAL*8 LATR,DECR
       REAL*8 COSENOH,COSENOA
       REAL*8 ALTR
       CHARACTER*1 ESTADO
C
       COMMON/BLKOBS/LAT,LONG
       COMMON/BLKOOC/HORTO,HOCASO,AORTO,AOCASO
C
       IF(LAT.GE.0)THEN
         IF(DECSUN.GE.90.-LAT)THEN
           ESTADO='C'
         ELSE IF(DECSUN.LE.LAT-90.)THEN
           ESTADO='I'
         ELSE
           ESTADO='O'
         END IF
       ELSE
         IF(DECSUN.LE.-(90.+LAT))THEN
           ESTADO='C'
         ELSE IF(DECSUN.GE.90.+LAT)THEN
           ESTADO='I'
         ELSE
           ESTADO='O'
         END IF
       END IF
C
       LATR=LAT*PI/180.D0
       DECR=DECSUN*PI/180.D0
       ALTR=(90.D0-DISZEN)*PI/180.D0
       COSENOH=(SIN(ALTR)-SIN(LATR)*SIN(DECR))/(COS(LATR)*COS(DECR))
       IF(COSENOH.LT.-1.) COSENOH=-1.
       IF(COSENOH.GT.+1.) COSENOH=+1.
       COSENOA=(SIN(ALTR)*SIN(LATR)-SIN(DECR))/(COS(ALTR)*COS(LATR))
       IF(COSENOA.LT.-1.) COSENOA=-1.
       IF(COSENOA.GT.+1.) COSENOA=+1.
       HOCASO=ACOS(COSENOH)
       HOCASO=HOCASO*180.D0/PI
       AOCASO=ACOS(COSENOA)
       AOCASO=AOCASO*180.D0/PI
       HORTO=360.D0-HOCASO
       AORTO=360.D0-AOCASO
       HORTO=HORTO/15.D0
       HOCASO=HOCASO/15.D0
C
       END
C **********************************************************************
C                                                     
       SUBROUTINE TIEMPOU(ARSUN,H,LONG,TS0,TSL,TU)
C
C Calculo del Tiempo Sidereo Local y del Tiempo Universal.
C
       IMPLICIT NONE
C
       REAL*8 ARSUN              ! Entrada
       REAL*8 H,LONG,TS0         ! Entrada
       REAL*8 TSL,TU             ! Salida
C
       REAL*8 DTS,TSG
C
       TSL=ARSUN+H
       IF(TSL.GT.24.)TSL=TSL-24.
       TSG=TSL-LONG
       IF(TSG.GT.24.)TSG=TSG-24.
       IF(TSG.LT.0.)TSG=TSG+24.
       DTS=TSG-TS0
       IF(DTS.LT.0.)DTS=DTS+24.
       TU=.997270*DTS
       IF(TU.LT.0.)TU=TU+24.
C
       END
C
C **********************************************************************
C     Pasamos un valor escrito como character formato +000 00 00 (II=1) 
C     o +00 00 00 (II=2) a decimal DD.dddd 
C                                                
       SUBROUTINE CH2DEC(CH11,DD,II)
C
       IMPLICIT NONE
       INTEGER II
       INTEGER A1,A2,A3
       REAL*8 DD
       CHARACTER*11 CH11
C
       IF(II.EQ.1) THEN
          IF(CH11(2:4).EQ.'  ') THEN
             A1=0
          ELSE
             READ(CH11(2:4),*) A1 ! Grados
          END IF
          IF(CH11(6:7).EQ.'  ') THEN
             A2=0
          ELSE
             READ(CH11(6:7),*) A2 ! Minutos
          END IF
          IF(CH11(9:10).EQ.'  ') THEN
             A3=0
          ELSE
             READ(CH11(9:10),*) A3 ! Segundos
          END IF
       ELSE IF (II.EQ.2) THEN
          IF(CH11(2:3).EQ.'  ') THEN
             A1=0
          ELSE
             READ(CH11(2:3),*) A1 ! Grados
          END IF
          IF(CH11(5:6).EQ.'  ') THEN
            A2=0
          ELSE
            READ(CH11(5:6),*) A2 ! Minutos
          END IF
          IF(CH11(8:9).EQ.'  ') THEN
             A3=0
          ELSE
             READ(CH11(8:9),*) A3 ! Segundos
          END IF
       END IF
       DD=DBLE(A1)+DBLE(A2)/60.D0+DBLE(A3)/3600.D0
       IF(CH11(1:1).EQ.'-') DD=-DD
       END
C
C **********************************************************************
C     Pasamos la hora de formato XX.xxxx (real) a  XXh XXm (character) 
C                                                        
      SUBROUTINE R2CH(RR,CH)
C
       INTEGER RR1,RR2
       REAL*8 RR
       CHARACTER*7 CH
C
       RR1=INT(RR)
       RR2=(RR-RR1)*60.
       WRITE(CH(1:2),'(I2)')RR1
       WRITE(CH(3:4),'(A1)')'h '
       WRITE(CH(5:6),'(I2)')RR2
       WRITE(CH(7:7),'(A1)')'m'
C
       END
C
C **********************************************************************
C     Pasamos coordenadas de formato XX.xxxx (real) a +gg mm ss (character) 
C                                                        
      SUBROUTINE RR2CH(RR,CH)
C
       INTEGER RR1,RR2,RR3
       REAL*8 RR
       CHARACTER*9 CH
C
       RR1=INT(RR)
       RR2=INT(ABS((RR-RR1)*60.))
       RR3=INT((ABS((RR-RR1)*60.)-RR2)*60.)
       WRITE(CH(1:3),'(I3)')RR1
       WRITE(CH(4:4),'(A1)') ' '
       WRITE(CH(5:6),'(I2)')RR2
       WRITE(CH(7:7),'(A1)')' '
       WRITE(CH(8:9),'(I2)')RR3
C
      END
C
C**********************************************************************
C     Pasamos de coordenadas ecuatoriales a horizontales
C
       SUBROUTINE ECU2HOR(TSL,LAT,AR,DEC,AZI,ALT)
C
       IMPLICIT NONE
       REAL*8 TSL,LAT           !Entrada
       REAL*8 AR,DEC          ! Entrada
       REAL*8 HOR
       REAL*8 AZI,ALT           ! Salida
       REAL*8 PI,G2R,R2G
       PARAMETER(PI=3.141592654D0)
       REAL*8 SENOAZ,COSEAZ
       REAL*8 HORR,DECR,LATR    ! Angulos en radianes
C
       G2R=PI/180.D0
       R2G=180.D0/PI
       HOR=TSL-AR    ! Angulo horario  (en horas)
       HORR=HOR*G2R*15.D0
       DECR=DEC*G2R
       LATR=LAT*G2R
       SENOAZ=COS(DECR)*SIN(HORR)
       COSEAZ=-(SIN(DECR)*COS(LATR))+COS(DECR)*SIN(LATR)*COS(HORR)
       AZI=ATAN2(SENOAZ,COSEAZ)
       AZI=AZI*R2G
       IF(AZI.LT.0.)AZI=AZI+360.D0
       ALT=ASIN(SIN(LATR)*SIN(DECR)+COS(LATR)*COS(DECR)*COS(HORR))
       ALT=ALT*R2G
C
      END
C
C**********************************************************************
C     Cambio de coordenadas ecuatoriales debido a la precesion de los
C     equinocios.

        SUBROUTINE CHANEQUI(TII,TFF,ARII,ARFF,DECII,DECFF)
C
	REAL*8  TII,TFF,TI,TF
	REAL*8 ARII,ARFF,DECII,DECFF
        REAL*8 X(3),X0(3)
        REAL*8 M(3,3)
        REAL*8 GIO,ZETA,TETA
	REAL*8 PI
	PARAMETER(PI=3.141592654D0)
C
        TI=(TII-2000.)/100.
        TF=(TFF-2000.-100.*TI)/100.

        GIO=((2306.2181D0+1.39656D0*TI-.000139D0*TI*TI)*TF+
     >  (.30188D0-.000344D0*TI)*TF*TF+.017998D0*TF*TF*TF)/3600.D0
        ZETA=GIO+((.7928+.00041*TI)*TF*TF+.000205*TF*TF*TF)/3600.D0
        TETA=((2004.3109D0-.8533D0*TI-.000217D0*TI*TI)*TF-
     >  (.42665D0+.000217D0*TI)*TF*TF-.041833*TF*TF*TF)/3600.D0
C
        GIO=GIO*PI/180.D0
        ZETA=ZETA*PI/180.D0
        TETA=TETA*PI/180.D0
C
        M(1,1)=-SIN(GIO)*SIN(ZETA)+COS(GIO)*COS(TETA)*COS(ZETA)
        M(1,2)=-COS(GIO)*SIN(ZETA)-SIN(GIO)*COS(ZETA)*COS(TETA)
        M(1,3)=-SIN(TETA)*COS(ZETA)
        M(2,1)=SIN(GIO)*COS(ZETA)+COS(GIO)*COS(TETA)*SIN(ZETA)
        M(2,2)=COS(GIO)*COS(ZETA)-SIN(GIO)*COS(TETA)*SIN(ZETA)
        M(2,3)=-SIN(TETA)*SIN(ZETA)
        M(3,1)=COS(GIO)*SIN(TETA)
        M(3,2)=-SIN(GIO)*SIN(TETA)
        M(3,3)=COS(TETA)
C
        X0(1)=COS(DECII*PI/180.D0)*COS(ARII*15.*PI/180.D0)
        X0(2)=COS(DECII*PI/180.D0)*SIN(ARII*15.*PI/180.D0)
        X0(3)=SIN(DECII*PI/180.D0)
C
        DO 505,K=1,3
          X(K)=0.
          DO 500,KK=1,3
            X(K)=X(K)+X0(KK)*M(K,KK)
500       CONTINUE
505     CONTINUE
        ARFF=ATAN2(X(2),X(1))
        DECFF=ASIN(X(3))
        ARFF=ARFF*180.D0/PI/15.D0
        IF(ARFF.LT.0.)ARFF=ARFF+24.D0
        DECFF=DECFF*180.D0/PI
	END
C
C**********************************************************************
C
      SUBROUTINE CULMINA(LONG,LAT,TS0,AR,DEC,ALTCUL,DIRCUL,TUCUL)
C
C     Se calcula la altura, direccion y tiempo universal de la culminacion
C     de una estrella a partir de sus AR. y DEC y la latitud del observador
C
      IMPLICIT NONE
C
      REAL*8 LONG,LAT
      REAL*8 AR,DEC
      REAL*8 ALTCUL,TUCUL
      REAL*8 TS0,TSL
C
      CHARACTER*1 DIRCUL
C
C     Comprobamos que el objeto sea visible
      IF(LAT.GE.0.) THEN
         IF(DEC.LE.LAT-90.) THEN
            DIRCUL='I'
            RETURN
         END IF
      ELSE
         IF(DEC.GE.90.+LAT) THEN
            DIRCUL='I'
            RETURN
         END IF
      END IF
C     Calculamos la altura y direccion de la culminacion
      IF(LAT.GE.0)THEN
        IF(DEC.LE.LAT)THEN
          ALTCUL=90.-LAT+DEC
          DIRCUL='S'
        ELSE
          ALTCUL=90.-DEC+LAT
          DIRCUL='N'
        END IF
      ELSE
        IF(DEC.LE.LAT)THEN
          ALTCUL=90.-ABS(DEC-LAT)
          DIRCUL='S'
        ELSE
          ALTCUL=90.-ABS(LAT-DEC)
          DIRCUL='N'
        END IF
      END IF
C     Calculamos T.U. del momento de la ocultacion
      CALL TIEMPOU(AR,0.D0,LONG,TS0,TSL,TUCUL)
C
      END
C
C**********************************************************************
C
       SUBROUTINE IOWIN1(IOTXT)
C
C     Escribimos el texto IOTXT en la ventanita de I/O y ayuda
C
       CHARACTER*40 IOTXT
C
       CALL PGSVP(0.0,1.0,0.0,1.0)
       CALL PGSWIN(0.0,2.0,0.0,1.0)
       CALL PGSLW(1)
       CALL PGSCI(16)
       CALL PGSFS(1)
       CALL PGRECT(0.60,1.50,0.00,0.06)
       CALL PGSCI(1)
       CALL PGSCF(3)
       CALL PGTEXT(0.61,0.02,IOTXT) 
       CALL PGSCF(2)
       CALL PGSFS(2)
C
       END
C**********************************************************************
C
       SUBROUTINE IOWIN2(IOTXT1,IOTXT2)
C
C     Escribimos el texto IOTXT en la ventanita de I/O y ayuda
C
       CHARACTER*50 IOTXT1,IOTXT2
C
       CALL PGSVP(0.0,1.0,0.0,1.0)
       CALL PGSWIN(0.0,2.0,0.0,1.0)
       CALL PGSLW(1)
       CALL PGSCI(16)
       CALL PGSFS(1)
       CALL PGRECT(0.60,1.50,0.00,0.06)
       CALL PGSCI(1)
       CALL PGSCF(3)
       CALL PGSCH(.8)
       CALL PGTEXT(0.62,0.035,IOTXT1) 
       CALL PGTEXT(0.62,0.010,IOTXT2) 
       CALL PGSCF(2)
       CALL PGSFS(2)
       CALL PGSCH(1.)
C
       END
C**********************************************************************
C
       SUBROUTINE MOONPOS(FJ,IYEAR,DAYS,MALFA,MDELTA,FASE)
C
C     Calculamos la posicion de la Luna siguiendo los pasos que se dan
C     'Practical Astronoy with your calculator' Peter Duffet-Smith,
C     cap.40 pg 91 y sgtes.
C     Tenemos como entrada un momento dado por la fecha juliana y como
C     salida las coord. ecuatoriales de la Luna en ese momento.
C     Elegimos como epoca de partida 1975.0 en la que :
C     Longitud de la Luna               MLONG   = 124.8756   grados
C     Longitud del perigeo              MPER    = 145.9601     "
C     Longitud del nodo                 MLONGN  = 248.6441     "
C     Inclinacion de orbita de Luna     MOINC   =   5.1453     "
C     Excentricidad orbita Luna         MOEXC   =   0.0549  
C     Tamanio angular medio de Luna     MSIZE   =   0.5181     "
C     Semieje mayor de orbita Luna      MSEJE   =   384401   Km
C     Paralaje a su distancia media     MPAR    =   0.9507   grados
C     Longitud ecliptica del Sol        SELONG  = 279.041470   "
C     Longitud ecliptica de su perigeo  SPELONG = 282.510396   "
C     Excentricidad de su orbita        SOEXC   =   0.016720
C
       IMPLICIT NONE
C
       INTEGER I
       INTEGER IYEAR
C
       REAL*8 MLONG,MPER,MLONGN,MOINC,MOEXC,MSIZE,MSEJE,MPAR
       REAL*8 DAYS,PI,G2R,R2G
       REAL*8 SELONG,SPELONG,SOEXC
       REAL*8 NN,EC
       REAL*8 SUNLONG  ! Longitud ecliptica del Sol en un momento dado
       REAL*8 SUNANOM  ! Anomalia media del Sol     "  "     "     "
       REAL*8 MOONLONG ! Longitud ecliptica de la Luna en un momento dado
       REAL*8 MOONLC,MOONLV  !     "        "     "  corregida Y verdadera
       REAL*8 MOONANOM ! Anomalia media     "  "   "   "       "     "
       REAL*8 MOONAC   ! Anomalia media Luna, corregida
       REAL*8 MOONLN   ! Longitud media del nodo ascendente
       REAL*8 MOONLNC  !   "        "        "      "       corregida
       REAL*8 EVEC     ! Evection
       REAL*8 ANECU    ! Ecuacion anual
       REAL*8 VARIA    ! Variacion
       REAL*8 ECC      ! Ecuacion del centro
       REAL*8 CCC
       REAL*8 MLANDA,MBETA ! Latitud y longitud eclipticas de la Luna
       REAL*8 FJ,MALFA,MDELTA  
       REAL*8 KK1,KK2
       REAL*8 FASE        !  Fraccion iluminada de cara Luna que vemos
       REAL*8 DANG        !  Angulo Tierra-Luna-Sol
C     
C     Inicializamos las variables
       PI=3.141592654D0
       G2R=PI/180.D0
       R2G=180.D0/PI
       MLONG   = 124.8756D0
       MPER    = 145.9601D0
       MLONGN  = 248.6441D0
       MOINC   =   5.1453D0
       MOEXC   =   0.0549D0
       MSIZE   =   0.5181D0
       MSEJE   =   384401D0
       MPAR    =   0.9507D0
       SELONG  = 279.041470D0
       SPELONG = 282.510396D0
       SOEXC   =   0.016720D0
C     Sumamos los dias de los anios desde 1975 hasta el actual.
       DO I=1975,IYEAR-1
          IF((I/4.-INT(I/4.)).EQ.0) THEN 
             DAYS=DAYS+366.D0
          ELSE
             DAYS=DAYS+365.D0
          END IF
       END DO
C
       NN=360.D0*DAYS/365.25D0
 30    IF(NN.GT.360.D0) THEN
          NN=NN-360.D0
          GOTO 30
       END IF
C     Calculamos la anomalia media del Sol
       SUNANOM=NN+SELONG-SPELONG
       IF(SUNANOM.LT.0.) SUNANOM=SUNANOM+360.D0
C     Calculamos la longitud ecliptica del Sol
       EC=(360.D0/PI)*SOEXC*DSIN(SUNANOM*PI/180.D0)
       SUNLONG= NN+EC+SELONG
       IF(SUNLONG.GT.360.D0) SUNLONG=SUNLONG-360.D0
C     Calculamos la longitud media de la Luna
       MOONLONG=360.D0*DAYS/27.3217D0+MLONG
C     Quitamos multiplos de 360, para que 0 > MOONLONG > 360 
 50    IF(MOONLONG.GT.360.D0) THEN
          MOONLONG=MOONLONG-360.D0
          GOTO 50
       ENDIF
C     Calculamos la anomalia media de la luna
       MOONANOM=MOONLONG-(360.D0/365.25D0)*(DAYS/8.85D0)-MPER
C     Aniadimos multiplos de 360, para que 0 > MOONANOM > 360 
 60    IF(MOONANOM.LT.0.) THEN
          MOONANOM=MOONANOM+360.D0
          GOTO 60
       ENDIF
C     Calculamos la longitud media del nodo ascendente
       MOONLN=MLONGN-(360.D0/365.25D0)*(DAYS/18.61D0)
C     Aniadimos multiplos de 360, para que 0 > MOONLN > 360 
       DO I=1,100
          IF(MOONLN.LT.0.) THEN
             MOONLN=MOONLN+360.D0
          ELSE
             GOTO 70
          ENDIF
       END DO
 70    CONTINUE
C     Calculamos la variacion de la excentricidad aparente de la orbita 
C     de la Luna (evection). Esta se le aplicara a la longitud y anomalia
C     medias de la Luna
       EVEC=1.274D0*DSIN((2.*(MOONLONG-SUNLONG)-MOONANOM)*G2R)
C     Calculamos la variacion de la distancia Sol-Tierra a lo largo del
C     anio (ecuacion anual). Esta se le restara a  la longitud y anomalia
C     medias de la Luna
       ANECU=0.18581D0*DSIN(SUNANOM*G2R)
C     Una tercera correccion que restar a la anomalia media de la Luna
       CCC=0.37D0*DSIN(SUNANOM*G2R)
C     Con estas tres correciones calculamos la anomalia media corregida
C     de la Luna
       MOONAC=MOONANOM+EVEC-ANECU-CCC
C     La ecuacion del centro sera:
       ECC=6.289D0*DSIN(MOONAC*G2R)
C     Calculamos la longitud corregida 
       MOONLC=MOONLONG+EVEC-ANECU+ECC
C     Calculamos lo que afecta el que la Luna se encuentre a un lado u otro
C     de la Tierra, es decir mas cerca o mas lejos del Sol, y por tanto
C     la variacion de la atraccion gravitatorio del Sol sobre la Luna
C     (variacion) Esta se le restara a la anomalia media de la Luna.
       VARIA=0.65831D0*DSIN(2.D0*(MOONLC-SUNLONG)*G2R)
C     Se la sumamos a la longitud corregida y nos da la longitud  verdadera
       MOONLV=MOONLC+VARIA
C     Con esta podemos calcular la longitud y latitud eclipticas.
       MOONLNC=MOONLN-0.16D0*DSIN(SUNANOM*G2R)
C     Para evitar ambiguidad en arcotangente, kk1 y kk2 deben ser muy parecidos
       KK2=MOONLV-MOONLNC
       KK1=DATAN(DTAN((KK2)*G2R)*DCOS(MOINC*G2R))
       KK1=KK1*R2G
 75    IF(ABS(KK2-KK1).GT.5.) THEN 
          IF(KK2.GT.KK1) THEN
             KK1=KK1+180.
          ELSE
             KK1=KK1-180.
          END IF
          GOTO 75
       END IF
C
       MLANDA=KK1+MOONLNC
       MBETA=DASIN(DSIN(KK2*G2R)*DSIN(MOINC*G2R))
       MBETA=MBETA*R2G
C     Pasamos a coord. ecuatoriales
       CALL ECLI2ECU(FJ,MLANDA,MBETA,MALFA,MDELTA)
C     Calculamos la fase de la Luna
       DANG=180.D0-SUNLONG+MOONLV
       FASE=0.5D0*(1.D0+COS(DANG*G2R))
C     
       END
c
C**********************************************************************
C
       SUBROUTINE NDAYS(IIYEAR,IIMONTH,IIDAY,IIHOUR,IIMIN,IISEC,DAYS)
C
C     Calculamos el numero de dias y fracc. (DAYS)  transcurridos del anio. 
C     Se considera Enero 0.0 a la medianoche entre 30 y 31 Dec del
C     anio anterior. Asi mediodia del dia 1 de Enero sera 1.5 
C
       IMPLICIT NONE
C
       INTEGER IIYEAR,IIMONTH,IIDAY,IIHOUR,IIMIN,IISEC
       INTEGER IMONTH,IDAYS
C
       REAL LEAP
       REAL*8 DAYS
C
       LOGICAL LLEAP    !TRUE -> Anio bisiesto
C
       LLEAP=.FALSE.
       LEAP=(IIYEAR/4.-INT(IIYEAR/4.))
       IF(LEAP.EQ.0.) LLEAP=.TRUE.
       IF(IIMONTH.GT.2) GOTO 80
       IMONTH=IIMONTH-1
       IF(LLEAP) THEN
          DAYS=IMONTH*62.D0
       ELSE
          DAYS=IMONTH*63.D0
       END IF
       DAYS=DAYS/2.D0
       IDAYS=INT(DAYS)
       GOTO 120
 80    IMONTH=IIMONTH+1
       DAYS=IMONTH*30.6D0
       IDAYS=INT(DAYS)
       IF(LLEAP) THEN
          IDAYS=IDAYS-62
       ELSE
          IDAYS=IDAYS-63
       END IF
 120   DAYS=IDAYS+IIDAY+IIHOUR/24.D0+IIMIN/1440.D0+IISEC/86400.D0
C
       END
C
C**********************************************************************
C     Pasamos de coordenadas eclipticas a ecuatoriales 
C
       SUBROUTINE ECLI2ECU(FJ,ECLONG,ECLAT,ALFA,DELTA)
C
       IMPLICIT NONE
C
       REAL*8 ECLONG,ECLAT  ! Latitud y longitud eclipticas (entrada)
       REAL*8 ALFA,DELTA    ! Coordenadas ecuatoriales A.R. Y DEC (salida) 
       REAL*8 FJ,EPSI
       REAL*8 G2R,R2G,PI
C
       PI=3.141592654D0
       G2R=PI/180.D0
       R2G=180.D0/PI
C     Calculamos la oblicuidad de la ecliptica (epsil nos la da en radianes)
       CALL EPSIL(FJ,EPSI)
C     Calculamos la declinacion
       DELTA=DASIN(DSIN(ECLAT*G2R)*DCOS(EPSI)
     >             +DCOS(ECLAT*G2R)*DSIN(EPSI)*DSIN(ECLONG*G2R))
       DELTA=DELTA*R2G
C     Calculamos la ascension recta
       ALFA=DATAN(DTAN(ECLONG*G2R)*DCOS(EPSI)
     >             -(DTAN(ECLAT*G2R)*DSIN(EPSI)/DCOS(ECLONG*G2R)))
       ALFA=ALFA*R2G
C     Para evitar ambiguedad al tomar el inverso de la tangente, nos 
C     aseguramos de que ALFA y ECLONG esten en el mismo cuadrante+
        IF(ECLONG.GT.90.)  ALFA=ALFA+180.
        IF(ECLONG.GT.270.) ALFA=ALFA+180.
C
       ALFA=ALFA/15.D0
C
       END
C
C**********************************************************************
C     Mostramos la ayuda correspondiente a la ventana que se haya elegido (NN)
C
      SUBROUTINE HELP(NN)
C
      IMPLICIT NONE
C
      INTEGER NN
C
      IF(NN.EQ.1) THEN
       CALL IOWIN2('Click left button to give the name of a new star  ',
     >             'or right button for a file with a list of objects ')
      ELSEIF(NN.EQ.2) THEN
       CALL IOWIN2('Enter right ascension for a new star              ',
     >             'with format :  00 00 00                           ')
      ELSEIF(NN.EQ.3) THEN
       CALL IOWIN2('Enter declination for a new star                  ',
     >             'with format : +00 00 00                           ')
      ELSEIF(NN.EQ.4) THEN
       CALL IOWIN2('Enter equinox for a new star                      ',
     >             'with format : 0000.00                             ')
      ELSEIF(NN.EQ.5) THEN
       CALL IOWIN2('Draw the position of the star with the coordinates',
     >             'shown for the chosen observatory                  ')
      ELSEIF(NN.EQ.6) THEN
       CALL IOWIN2('Delete the traces of the stars and the Moon       ',
     >             'currently drawed                                  ')
      ELSEIF(NN.EQ.7) THEN
       CALL IOWIN2('Click left mouse button to change the name or the ',
     >             'right one to select another observatory           ')
      ELSEIF(NN.EQ.8) THEN
       CALL IOWIN2('Enter the longitude of a new observatory          ',
     >             'with format : +000 00 00                          ')
      ELSEIF(NN.EQ.9) THEN
       CALL IOWIN2('Enter the latitude of a new observatory           ',
     >             'with format : +00 00 00                           ')
      ELSEIF(NN.EQ.10) THEN
       CALL IOWIN2('Enter the height of a new observatory             ',
     >             'with format : 0000.                               ')
      ELSEIF(NN.EQ.11) THEN
       CALL IOWIN2('Enter a new date                                  ',
     >             'with format : DD-MM-YYYY                          ')
      ELSEIF(NN.EQ.12) THEN
       CALL IOWIN2('Create a file "ves.lis" that list the drawed stars',
     >             'This can be used later as input-file              ')
      ELSEIF(NN.EQ.13) THEN
       CALL IOWIN2('Create a postscript file with the current plot and',
     >             'data like coord., culmination, visibility, etc ...')
      ELSEIF(NN.EQ.14) THEN
       CALL IOWIN2('Yes, exactly. Exit = Quit = Bye                   ',
     >             'Adios, auf Wiedersehen ....                       ')
      ELSEIF(NN.EQ.15) THEN
       CALL IOWIN1(' Set Elevation for Y-axis                ')
      ELSEIF(NN.EQ.16) THEN
       CALL IOWIN1(' Set sec-Z ~ airmass for Y-axis          ')
      ELSEIF(NN.EQ.17) THEN
       CALL IOWIN1(' Set U.T. for X-axis                     ')
      ELSEIF(NN.EQ.18) THEN
       CALL IOWIN1(' Plot elevation or sec-Z vs azimut       ')
      END IF
C
      END
