Changeset 672 for trunk/LMDZ.MARS/libf

Timestamp:

May 25, 2012, 3:43:44 PM (13 years ago)

Author:

tnavarro

Message:

change perihelion date in solar longitude for any orbit, not only present day one

File:

• trunk/LMDZ.MARS/libf/phymars/tabfi.F (modified) (3 diffs)

trunk/LMDZ.MARS/libf/phymars/tabfi.F

@@ -445 +445 @@
           if(ierr.ne.0) goto 1160
           write(*,*) 
-          write(*,*) ' peri_ls asked :',peri_ls
-          call lsp2solp(peri_ls,peri_day)
-          write(*,*) ' peri_day (new value):',peri_day
+          write(*,*) 'peri_ls asked:',peri_ls
+          write(*,*) 'for aphelion =',aphelie
+          write(*,*) 'perihelion =',periheli
+          write(*,*) 'and',year_day,'sols/year'
+          call lsp2solp(peri_ls,peri_day,aphelie,periheli,year_day)
+          write(*,*) 'peri_day (new value):',peri_day
 
 
@@ -549 +552 @@
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 ! gives sol at perihelion for ls at perihelion (for precession cycles)
-      subroutine lsp2solp(lsp,solp)
+      subroutine lsp2solp(lsp,solp,aphelie,periheli,year_day)
 
       implicit none
@@ -555 +558 @@
       real lsp     ! Input: ls at perihelion
       real solp    ! Output: sol at perihelion
-      
+      real aphelie,periheli,year_day ! Input: parameters
+ 
 !  Local:
       double precision zx0 ! eccentric anomaly at Ls=0
-      double precision year_day 
       double precision e_elips
       double precision pi,degrad 
       
-      parameter (year_day=668.6d0) ! number of sols in a martian year
-      parameter (e_elips=0.0934d0)  ! eccentricity of orbit
       parameter (pi=3.14159265358979d0)
       parameter (degrad=57.2957795130823d0)
-      
+
+      e_elips=(aphelie-periheli)/(aphelie+periheli)      
       zx0 = -2.0*datan(dtan(0.5*lsp/degrad)
      .          *dsqrt((1.-e_elips)/(1.+e_elips)))