source: trunk/LMDZ.GENERIC/libf/phystd/stellarlong.F @ 603

Last change on this file since 603 was 590, checked in by aslmd, 13 years ago

LMDZ.GENERIC: Introduced global1d in callcorrk so that global (using sza) or local (using latitude) 1D simulations can be carried out. Converted all astronomical distances in AU instead of Mkm. This might cause problems with old start files. So added a test in iniorbit. A quite dirty test, but thatll do the job.

File size: 2.4 KB
Line 
1      SUBROUTINE stellarlong(pday,pstellong)
2      IMPLICIT NONE
3
4c=======================================================================
5c
6c   Objet:
7c   ------
8c
9c      Calcul de la distance soleil-planete et de la declinaison
10c   en fonction du jour de l'annee.
11c
12c
13c   Methode:
14c   --------
15c
16c      Calcul complet de l'elipse
17c
18c   Interface:
19c   ----------
20c
21c      Uncommon comprenant les parametres orbitaux.
22c
23c   Arguments:
24c   ----------
25c
26c   Input:
27c   ------
28c   pday          jour de l'annee (le jour 0 correspondant a l'equinoxe)
29c
30c   Output:
31c   -------
32c   pdist_star     distance entre le soleil et la planete
33c                 ( en unite astronomique pour utiliser la constante
34c                  solaire terrestre 1370 Wm-2 )
35c   pdecli        declinaison ( en radians )
36c
37c=======================================================================
38c-----------------------------------------------------------------------
39c   Declarations:
40c   -------------
41
42#include "planete.h"
43#include "comcstfi.h"
44
45c arguments:
46c ----------
47
48      REAL pday,pdist_star,pdecli,pstellong
49      LOGICAL lwrite
50
51c Local:
52c ------
53
54      REAL zanom,xref,zx0,zdx,zteta,zz
55      INTEGER iter
56
57
58c-----------------------------------------------------------------------
59c calcul de l'angle polaire et de la distance au soleil :
60c -------------------------------------------------------
61
62c  calcul de l'zanomalie moyenne
63
64      zz=(pday-peri_day)/year_day
65      pi=2.*asin(1.)
66      zanom=2.*pi*(zz-nint(zz))
67      xref=abs(zanom)
68
69c  resolution de l'equation horaire  zx0 - e * sin (zx0) = xref
70c  methode de Newton
71
72      zx0=xref+e_elips*sin(xref)
73      DO 110 iter=1,10
74         zdx=-(zx0-e_elips*sin(zx0)-xref)/(1.-e_elips*cos(zx0))
75         if(abs(zdx).le.(1.e-7)) goto 120
76         zx0=zx0+zdx
77110   continue
78120   continue
79      zx0=zx0+zdx
80      if(zanom.lt.0.) zx0=-zx0
81
82c zteta est la longitude solaire
83
84      zteta=2.*atan(sqrt((1.+e_elips)/(1.-e_elips))*tan(zx0/2.))
85
86      pstellong=zteta-timeperi
87
88      IF(pstellong.LT.0.) pstellong=pstellong+2.*pi
89      IF(pstellong.GT.2.*pi) pstellong=pstellong-2.*pi
90c-----------------------------------------------------------------------
91c   sorties eventuelles:
92c   ---------------------
93
94c     IF (lwrite) THEN
95c        PRINT*,'jour de l"annee   :',pday
96c        PRINT*,'distance au soleil (en unite astronomique) :',pdist_star
97c        PRINT*,'declinaison (en degres) :',pdecli*180./pi
98c     ENDIF
99
100      RETURN
101      END
Note: See TracBrowser for help on using the repository browser.