* 
*     PROGRAMME CALDAT    GF9901 
* 
* --- Objet ------------------------------------------------------------ 
* 
*     Consultation des equivalences entre systemes de datation. 
* 
*     Date du calendrier (jour, mois, an) 
*     Date Julienne Modifiee MJD 
*     Date Julienne - 2000000.5 
*     Annee besselienne. 
*     Annee et Numero du jour. 
* 
* --- Parametres ------------------------------------------------------- 
* 
*     an1  : annee limite inferieure (entier). 
*     an2  : annee limite superieure (entier). 
*     ab0  : annee origine des dates besseliennes (reel). 
*     djb0 : date julienne origine des dates besseliennes (reel). 
* 
*     L'annee an1 doit toujours etre inferieure a an2. 
*     Les dates juliennes modifiees MJD ont pour origine 2400000.5. 
*     La date besselienne ab0 doit etre coherente avec djb0. 
* 
*     Les valeurs ab0=1900.0 et djb0=2415020.31352d0 sont celles 
*     couramment utilisees pour calculer la date julienne correspondant 
*     au debut d'une annee besselienne et inversement. 
* 
* --- Remarque --------------------------------------------------------- 
* 
*     Quand l'utilisateur entre une date du calendrier (jour, mois, an) 
*     ou une annee avec le numero du jour, ces dates sont considerees 
*     a 0h pour calculer les dates juliennes et les annees besseliennes 
*     correspondantes. 
* 
*     Cependant, quand l'utilisateur entre une date julienne modifiee 
*     ou une annee besselienne, il n'est precise ni l'heure, ni la 
*     fraction de jour correspondantes dans les autres datations. 
* 
*     Les dates juliennes modifiees sont enregistrees et affichees avec 
*     2 decimales apres la virgule -> precision de 15 minutes environ. 
*     Les annees besseliennes sont enregistrees et affichees avec 
*     5 decimales apres la virgule -> precision de 5 minutes environ. 
* 
*     Cela signifie que si l'utilisateur entre une date julienne 
*     modifiee, la date besselienne correspondante peut avoir une 
*     incertitude d'environ 0.00003, alors que s'il entre une annee 
*     besselienne, les dates juliennes correspondantes sont justes 
*     a 0.01 pres. 
* 
* --- Declarations ----------------------------------------------------- 
* 
      implicit double precision (a-h,o-z) 
      integer an,an1,an2 
      character rep 
      character*4 mois(12) 
      character*8 nomj(7) 
* 
      dimension knom(7) 
* 
      parameter (an1=1750,an2=2250) 
      parameter (djn0=2000000.5d0) 
      parameter (ab0=1900.0,djb0=2415020.31352d0) 
* 
* --- Initialisation --------------------------------------------------- 
* 
      data antr/365.242198781d0/ 
* 
      data mois/ 
     .'JANV','FEVR','MARS','AVR.','MAI ','JUIN', 
     .'JUIL','AOUT','SEPT','OCT.','NOV.','DEC.'/ 
* 
      data nomj/ 
     .'   LUNDI','   MARDI','MERCREDI','   JEUDI', 
     .'VENDREDI','  SAMEDI','DIMANCHE'/ 
* 
      data knom /4, 4, 1, 4, 1, 3, 1/ 
* 
      call DATEDJ (1, 1, an1, 0, 0, 0, dj1) 
      djm1 = dj1 - 2400000.5d0 
      djn1 = dj1 - djn0 
      ab1 = an1 
      nj1 = 1 
* 
      call DATEDJ (31, 12, an2, 0, 0, 0, dj2) 
      djm2 = dj2-2400000.5d0 
      djn2 = dj2-djn0 
      ab2 = an2 
      call DATEDJ (1, 1, an2, 0, 0, 0, dj0) 
      nj2=dj2-dj0+1 
* 
* --- Faire un choix de datation --------------------------------------- 
* 
100   continue 
* 
      dj=0.d0 
* 
      call EFFECR 
* 
      write (*,1000) an1,an2,djm1,djm2,djn1,djn2,ab1,ab2,an1,nj1,an2,nj2 
      read (*,'(a)') rep 
* 
* --- Lire une date dans la datation choisie --------------------------- 
* 
      select case (rep) 
* 
      case ('0') 
         stop 
* 
      case ('1') 
         write (*,1001) 
         read (*,*,iostat=nerr) j,m,an 
         if (nerr.eq.0) then 
            if (an.ge.0.and.an.le.99) then 
               if (an.ge.50) an=1900+an 
               if (an.lt.50) an=2000+an 
               call DATEDJ (j,m,an,0,0,0,dj) 
               if (dj.lt.1.d-8) nerr=1 
            else 
               if (an.ge.an1.and.an.le.an2) then 
                  call DATEDJ (j,m,an,0,0,0,dj) 
                  if (dj.lt.1.d-8) nerr=1 
               else 
                  nerr=1 
               endif 
            endif 
         else 
            nerr=1 
         endif 
* 
      case ('2') 
         write (*,1002) 
         read (*,*,iostat=nerr) djm 
         if (nerr.eq.0.and.djm.ge.djm1.and.djm.le.djm2) then 
            d=dint(abs(djm)*1.d2+0.5d0)/1.d2 
            if (djm.lt.0.d0) d=-d 
            dj=d+2400000.5d0 
         else 
            nerr=1 
         endif 
* 
      case ('3') 
         write (*,1003) 
         read (*,*,iostat=nerr) djn 
         if (nerr.eq.0.and.djn.ge.djn1.and.djn.le.djn2) then 
            d=dint(abs(djn)*1.d2+0.5d0)/1.d2 
            if (djn.lt.0.d0) d=-d 
            dj=d+djn0 
         else 
            nerr=1 
         endif 
* 
      case ('4') 
         write (*,1004) 
         read (*,*,iostat=nerr) ab 
         if (nerr.eq.0.and.ab.ge.ab1.and.ab.le.ab2) then 
            ab=dint(ab*1.d5+0.5d0)/1.d5 
            dj=djb0+(ab-ab0)*antr 
         else 
            nerr=1 
         endif 
* 
      case ('5') 
         write (*,1005) 
         read (*,*,iostat=nerr) an,nj 
         if (nerr.eq.0) then 
            if (an.ge.0.and.an.le.99.and. 
     .          nj.ge.0.and.nj.le.366) then 
               if (an.ge.50) an=1900+an 
               if (an.lt.50) an=2000+an 
               call DATEDJ (1,1,an,0,0,0,dj) 
               if (dj.lt.1.d-8) then 
                  nerr=1 
               else 
                  dj=dj+nj-1 
               endif 
            else 
               if (an.ge.an1.and.an.le.an2.and. 
     .             nj.ge.0.and.nj.le.366) then 
                  call DATEDJ (1,1,an,0,0,0,dj) 
                  if (dj.lt.1.d-8) then 
                     nerr=1 
                  else 
                     dj=dj+nj-1 
                  endif 
               else 
                  nerr=1 
               endif 
            endif 
         else 
            nerr=1 
         endif 
* 
      case default 
         nerr=1 
* 
      endselect 
* 
* --- Erreur de choix -------------------------------------------------- 
* 
      if (nerr.ne.0) then 
         write (*,1006) 
         read (*,*) 
         goto 100 
      endif 
* 
* --- Calcul de la date dans les differentes datation ------------------ 
* 
      call DJDATE (dj,0,j,m,an,nh,nm,ns) 
      djm=dj-2400000.5d0 
      djn=dj-djn0 
      ab=ab0+(dj-djb0)/antr 
      call DATEDJ (1,1,an,0,0,0,djz) 
      nj=dj-djz+1 
      kj=JSEM(dj) 
      jj=knom(kj) 
* 
      write (*,1007) nomj(kj)(jj:8),j,mois(m),an,djm,djn,ab,an,nj,dj 
* 
* --- Encore ou Stop ? ------------------------------------------------- 
* 
      write (*,1008) 
      read (*,'(a)') rep 
      if (rep.ne.'N'.and.rep.ne.'n') goto 100 
      stop 
* 
* --- Formats ---------------------------------------------------------- 
* 
1000  format (/2x,'DATEUR : EQUIVALENCE ENTRE SYSTEMES DE DATATON'/ 
     .        /2x,'1 : Jour Mois An               ', 
     .         4x,'(entre 01 01 ',i4,' et 31 12 ',i4,')' 
     .        /2x,'2 : Date julienne modifiee MJD ', 
     .         4x,'(entre ',f10.2,' et ',f10.2,')', 
     .        /2x,'3 : Date julienne - 200000.5   ', 
     .         4x,'(entre ',f10.2,' et ',f10.2,')', 
     .        /2x,'4 : Annee Besselienne          ', 
     .         4x,'(entre ',f10.5,' et ',f10.5,')', 
     .        /2x,'5 : An et Numero du jour       ', 
     .         4x,'(entre ',i4,'  J',i3.3,' et ',i4,'  J',i3.3,')' 
     .        /2x,'0 : Sortie'/ 
     .        /2x,'Entrez votre choix : '$) 
1001  format (/2x,'Entrez une date calendaire (jour mois an) : ',$) 
1002  format (/2x,'Entrez une date julienne modifiee MJD : ',$) 
1003  format (/2x,'Entrez une date julienne - 2000000.5 : ',$) 
1004  format (/2x,'Entrez une annee besselienne : ',$) 
1005  format (/2x,'Entrez une annee et le numero du jour : ',$) 
1006  format (/2x,'Erreur de choix : Appuyez sur Entree ') 
1007  format (/2x,a,1x,i2.2,1x,a,1x,i4/ 
     .        /2x,'Date julienne modifiee MJD : ',f10.2 
     .        /2x,'Date julienne - 2000000.5  : ',f10.2 
     .        /2x,'Annee besselienne          : ',f10.5 
     .        /2x,'An et Numero du jour       : ',i4,'  J',i3.3 
     .        /2x,'Date Julienne              : ',f10.2) 
1008  format (/2x,'Voulez-vous faire un autre calcul (o/n) ? ',$) 
* 
      end 
* 
* 
* 
      subroutine DJDATE (tdj,icode,jour,mois,kan,nh,nm,ns) 
*----------------------------------------------------------------------- 
* 
*     Ref : GF9803 
* 
* --- Objet ------------------------------------------------------------ 
* 
*     Changement de calendrier : Date julienne -> Gregorien/Julien. 
* 
* --- Entree ------------------------------------------------------------ 
* 
*     tdj         Date julienne (reel dp). 
* 
*     icode       Indice code de la conversion (entier) 
*                 icode=0 conversion tdj -> jour-mois-kan et 00-00-00. 
*                 icode=1 conversion tdj -> jour-mois-kan et nh-00-00. 
*                 icode=2 conversion tdj -> jour-mois-kan et nh-nm-00. 
*                 icode=3 conversion tdj -> jour-mois-kan et nh-nm-ns. 
* 
* --- Sortie ----------------------------------------------------------- 
* 
*     jour        Quantieme du mois (entier). 
*     mois        Numero du mois (entier). 
*     kan         Millesime astronomique de l'annee (entier). 
* 
*     nh          Nombre d'heures (entier). 
*     nm          Nombre de minutes (entier). 
*     ns          Nombre de secondes (entier). 
* 
* --- Remarque --------------------------------------------------------- 
* 
*     Date du calendier (Gregorien/Julien) : 
*     Le calendier Gregorien commence le 15/10/1582 a 0h (DJ2299160.5). 
*     Avant cette date on utilise le calendrier Julien. 
* 
*     On utilise le formulaire de J. Meeus. 
*     La date julienne tdj doit etre superieure ou egale a 0. 
*     Si tdj<0, alors jour=mois=kan=nh=nm=ns=0. 
* 
*     L'annee kan est comptee en millesimes astronomiques. 
*     L'annee 1 avant J.C. des historiens correspond a l'an 0 des 
*     astronomes et ainsi de suite. 
*     Ex : l'an 80 avant J.C. correspond a l'an astronomique -79. 
* 
* --- Declarations ----------------------------------------------------- 
* 
      implicit double precision (a-h,o-z) 
* 
* --- Initialisation --------------------------------------------------- 
* 
      jour=0 
      mois=0 
      kan=0 
      nh=0 
      nm=0 
      ns=0 
* 
      if (tdj.lt.0.d0.or.icode.lt.0.or.icode.gt.3) return 
      if (icode.eq.0) dt=0.0d0 
      if (icode.eq.1) dt=0.5d0/24.d0 
      if (icode.eq.2) dt=0.5d0/1440.d0 
      if (icode.eq.3) dt=0.5d0/86400.d0 
* 
* --- Conversion ------------------------------------------------------- 
* 
      t=tdj+dt+0.5d0 
      z=int(t) 
      f=t-z 
* 
      if (z.lt.2299161.d0) then 
         a=z 
      else 
         alpha=int((z-1867216.25d0)/36524.25d0) 
         a=z+1.d0+alpha-int(alpha/4.d0) 
      endif 
* 
      b=a+1524.d0 
      c=int((b-122.1d0)/365.25d0) 
      d=int(365.25d0*c) 
      e=int((b-d)/30.6001d0) 
* 
      jour=b-d-int(30.6001d0*e) 
* 
      if (e.lt.13.5d0) then 
         mois=e-1.d0 
      else 
         mois=e-13.d0 
      endif 
* 
      if (mois.lt.2.5d0) then 
         kan=c-4715.d0 
      else 
         kan=c-4716.d0 
      endif 
* 
      f=f*24.d0 
      nh=f 
      f=(f-nh)*60.d0 
      nm=f 
      f=(f-nm)*60.d0 
      ns=f 
* 
* --- Sortie ----------------------------------------------------------- 
* 
      return 
      end 
* 
* 
* 
      subroutine DATEDJ (jour,mois,kan,nh,nm,ns,tdj) 
*----------------------------------------------------------------------- 
* 
*     Ref : GF9803 
* 
* --- Objet ------------------------------------------------------------ 
* 
*     Changement de calendrier : Gregorien/Julien -> Date julienne. 
* 
* --- Entree ------------------------------------------------------------ 
* 
*     jour        Quantieme du mois (entier). 
*     mois        Numero du mois (entier). 
*     kan         Millesime astronomique de l'annee (entier). 
* 
*     nh          Nombre d'heures (entier). 
*     nm          Nombre de minutes (entier). 
*     ns          Nombre de secondes (entier). 
* 
* --- Sortie ----------------------------------------------------------- 
* 
*     tdj         Date julienne (reel dp). 
* 
* --- Remarque --------------------------------------------------------- 
* 
*     Date du calendier (Gregorien/Julien) : 
*     Le calendier Gregorien commence le 15/10/1582 a 0h (DJ2299160.5). 
*     Avant cette date on utilise le calendrier Julien. 
* 
*     On utilise le formulaire de J. Meeus. 
*     Cette fonction n'est valable que pour les dates posterieures au 
*     1/1/4713 12h avant J.C. 
* 
*     L'annee kan est comptee en millesimes astronomiques. 
*     L'annee 1 avant J.C. des historiens correspond a l'an 0 des 
*     astronomes et ainsi de suite. 
*     Ex : l'an 80 avant J.C. correspond a l'an astronomique -79. 
* 
* --- Declarations ----------------------------------------------------- 
* 
      implicit double precision (a-h,o-z) 
* 
* --- Initialisation --------------------------------------------------- 
* 
      tdj=0.d0 
* 
      if ((kan.lt.-4713) 
     .   .or.(mois.lt.1.or.mois.gt.12) 
     .   .or.(jour.lt.1.or.jour.gt.31) 
     .   .or.(nh.lt.0.or.nh.gt.24) 
     .   .or.(nm.lt.0.or.nm.gt.59) 
     .   .or.(ns.lt.0.or.ns.gt.59)) return 
* 
      a=0.d0 
      b=0.d0 
      c=0.d0 
* 
* --- Conversion ------------------------------------------------------- 
* 
      if (mois.gt.2) then 
         y=kan 
         m=mois 
      else 
         y=kan-1.d0 
         m=mois+12 
      endif 
* 
      if (y.lt.0) then 
         c=-0.75d0 
      else 
         idate=kan*10000+mois*100+jour 
         if (idate.ge.15821015) then 
            a=int(y/100.d0) 
            b=2.d0-a+int(a/4.d0) 
         endif 
      endif 
* 
      tdj=int(365.25d0*y+c) 
     .   +int(30.6001d0*(m+1)) 
     .   +jour+b 
     .   +dfloat(nh)/24.d0 
     .   +dfloat(nm)/1440.d0 
     .   +dfloat(ns)/86400.d0 
     .   +1720994.5d0 
* 
* --- Sortie ----------------------------------------------------------- 
* 
      return 
      end 
* 
* 
* 
      function JSEM (tdj) 
*----------------------------------------------------------------------- 
* 
*     Ref : GF9803 
* 
* --- Objet ------------------------------------------------------------ 
* 
*     Calcul du jour de la semaine. 
* 
* --- Entree ------------------------------------------------------------ 
* 
*     tdj         Date Julienne (reel dp). 
* 
* --- Sortie ----------------------------------------------------------- 
* 
*     jsem        Jour de la semaine (entier). 
*                 jsem=1 lundi, jsem=2 mardi, ... , jsem=7 dimanche. 
* 
* --- Remarque --------------------------------------------------------- 
* 
*     Une date julienne se terminant par 0.5 est consideree comme le 
*     jour suivant a 0h. 
* 
*----- Calcul ---------------------------------------------------------- 
* 
      implicit double precision (a-h,o-z) 
* 
      JSEM=0 
* 
      if (tdj.lt.0.0d0) return 
* 
      k=tdj+0.5d0/86400.d0+1.5d0 
      k=mod(k,7) 
      if (k.eq.0) k=7 
* 
      JSEM=k 
* 
      return 
      end 
* 
* 
* 
      subroutine EFFECR 
*----------------------------------------------------------------------- 
* 
*     Ref : GF9802 
* 
* --- Objet ------------------------------------------------------------ 
* 
*     Effacement de l'ecran. 
* 
      call SYSTEM ('clear') 
* 
      return 
      end 
      