MODULE recomp_orb_param_mod
!=======================================================================
!
! Purpose: Recompute orbit parameters based on values by Laskar et al.,
! 2004
!
! Authors: RV 
!=======================================================================
      IMPLICIT NONE

      CONTAINS

      SUBROUTINE recomp_orb_param(final_iter)

      USE temps_mod_evol, ONLY: year_bp_ini, year_PEM, var_obl, var_ex, var_lsp
#ifndef CPP_STD      
      USE comconst_mod, ONLY: pi
      USE planete_h, ONLY: e_elips, obliquit, timeperi
#else
      use planete_mod, only: e_elips, obliquit, timeperi
      USE comcstfi_mod, only: pi

#endif

      USE lask_param_mod, only: yearlask,oblask,exlask,lsplask, &
   				end_lask_param_mod,last_ilask

      IMPLICIT NONE

!--------------------------------------------------------
! Input Variables
!--------------------------------------------------------

      integer,intent(in) :: final_iter ! Number of year iteration of the PEM 

!--------------------------------------------------------
! Output Variables
!--------------------------------------------------------

!--------------------------------------------------------
! Local variables 
!--------------------------------------------------------

      real :: Year	                ! Year of the simulation
      real :: timeperi_ls	        ! time of peri in ls
      integer ilask			! Loop variable


      ! **********************************************************************
      ! 0. Initializations
      ! **********************************************************************


          Year=year_bp_ini+year_PEM+final_iter

          timeperi_ls=timeperi*360/2/pi

!          Year=-953200
          print *, "recomp_orb_param, Old year=", year_bp_ini+year_PEM
          print *, "recomp_orb_param, New year=", Year
          print *, "recomp_orb_param, Old obl=", obliquit
          print *, "recomp_orb_param, Old ex=", e_elips
          print *, "recomp_orb_param, Old lsp=", timeperi_ls
          print *, "recomp_orb_param, Old timeperi=", timeperi

        do ilask=last_ilask,1,-1
           if(yearlask(ilask) .GT.Year) then
             if(var_obl) then
               obliquit=oblask(ilask+1)+(oblask(ilask)-oblask(ilask+1))*(Year-yearlask(ilask+1))/(yearlask(ilask)-yearlask(ilask+1))
             endif
             if(var_ex) then
               e_elips=exlask(ilask+1)+(exlask(ilask)-exlask(ilask+1))*(Year-yearlask(ilask+1))/(yearlask(ilask)-yearlask(ilask+1))
             endif
             if(var_lsp) then
               timeperi_ls=lsplask(ilask+1)+(lsplask(ilask)-lsplask(ilask+1))*(Year-yearlask(ilask+1))/(yearlask(ilask)-yearlask(ilask+1))
             endif
              exit
           endif
        enddo
        timeperi=timeperi_ls*2*pi/360

       print *, "recomp_orb_param, Final year of the PEM run:", Year
       print *, "recomp_orb_param, New obl=", obliquit
       print *, "recomp_orb_param, New ex=", e_elips
       print *, "recomp_orb_param, New timeperi=", timeperi

      END SUBROUTINE recomp_orb_param

!********************************************************************************   
      
      END MODULE recomp_orb_param_mod