source: trunk/LMDZ.GENERIC/libf/phystd/orbite.F @ 2613

Last change on this file since 2613 was 1384, checked in by emillour, 10 years ago

Generic GCM:

  • Some code cleanup: turning comcstfi.h into module comcstfi_mod.F90

EM

File size: 1.8 KB
Line 
1      subroutine orbite(pls,pdist_star,pdecli,pright_ascenc)
2
3      use planete_mod, only: p_elips, e_elips, timeperi, obliquit
4      use comcstfi_mod, only: pi
5      implicit none
6!==================================================================
7!     
8!     Purpose
9!     -------
10!     Distance from star and declination as a function of the stellar
11!     longitude Ls
12!     
13!     Inputs
14!     ------
15!     pls          Ls
16!
17!     Outputs
18!     -------
19!     pdist_star    Distance Star-Planet in UA
20!     pdecli        declinaison ( in radians )
21!     pright_ascenc right ascension ( in radians )
22!
23!=======================================================================
24
25c   Declarations:
26c   -------------
27
28c arguments:
29c ----------
30
31      REAL pday,pdist_star,pdecli,pright_ascenc,pls,i
32
33c-----------------------------------------------------------------------
34
35c Star-Planet Distance
36
37      pdist_star = p_elips/(1.+e_elips*cos(pls+timeperi))
38
39c Stellar declination
40
41c ********************* version before 01/01/2000 *******
42
43      pdecli = asin (sin(pls)*sin(obliquit*pi/180.))
44
45c********************* version after 01/01/2000 *******
46c     i=obliquit*pi/180.
47c     pdecli=asin(sin(pls)*sin(i)/sqrt(sin(pls)**2+
48c    & cos(pls)**2*cos(i)**2))
49c ******************************************************
50
51c right ascencion
52      If((pls.lt.pi/2.d0)) then
53         pright_ascenc= atan(tan(pls)*cos(obliquit*pi/180.))
54      else if((pls.gt.pi/2.d0).and.(pls.lt.3.d0*pi/2.d0)) then
55         pright_ascenc= pi+atan(tan(pls)*cos(obliquit*pi/180.))
56      else if((pls.gt.3.d0*pi/2.d0)) then
57         pright_ascenc= 2.d0*pi+atan(tan(pls)*cos(obliquit*pi/180.))
58      else if (Abs(pls-pi/2.d0).le.1.d-10) then
59         pright_ascenc= pi/2.d0
60      else
61         pright_ascenc=-pi/2.d0
62      end if
63         
64      RETURN
65      END
Note: See TracBrowser for help on using the repository browser.