source: trunk/SOURCES/Greeneem_files/lect-clim-act-greeneem_mois_lapsecouche_mod.f90 @ 237

module lect_clim_act_greeneem_mois_lapsecouche                    


use module3D_phy                        
use printtable
use declare_month_lapsecouche
use climat_perturb_mois_lapsecouche_mod


implicit none



contains


subroutine input_climat_ref()           ! routine qui permet d'initialiser les variables climatiques




real,dimension(nx,ny,mois) :: Tm       ! temperature mensuelle sur la calotte de ref
real,dimension(nx,ny,mois) :: Pm       ! precipitation mensuelle sur la calotte de ref


character(len=200)  :: filin                            ! nom temporaire
character (len=200) :: filtr_topo, filtr_t, filtr_p
real :: T_surf_ref                                      ! variable de travail
real :: zel                                             ! variable de travail
real :: ti, tj                     ! variable provisoire de lecture


namelist/clim_init_mois/filtr_topo,filtr_t,filtr_p    

! lecture par namelist 
!---------------------


! formats pour les ecritures dans 42
428 format(A)

rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
read(num_param,clim_init_mois)


write(num_rep_42,428)'!___________________________________________________________' 
write(num_rep_42,428) '&clim_perturb_mois                            ! module climat_perturb_mois_mod'
write(num_rep_42,'(A,A)')   'filtr_topo =', filtr_topo
write(num_rep_42,'(A,A)')   'filtr_t =', filtr_t
write(num_rep_42,'(A,A)')   'filtr_p =', filtr_p
write(num_rep_42,*)'/'                      
write(num_rep_42,428) '! fichiers temperature et precip : 12 mois'
write(num_rep_42,428) '! faire un bloc namelist par snapshot'
write(num_rep_42,*)



call lect_lapserate_months       ! lit les lasperate
                                 ! version actuelle : on utilise un
                                 ! lapserate parametre homogene annuel
                                 ! on le passe en mensuel avec un cosinus



! lecture des fichiers snapshots pour tout geoplace
! -------------------------------------------------
write(6,*) 'fichiers snapshots'

    ! lecture de la topo de surface de reference du GCM
   filin=trim(dirnameinp)//trim(filtr_topo)
   write(6,*) filin
   open(20,file=filin)
   do j=1,ny 
      do i=1,nx
         read (20,*) S0CLIM(i,j)
         S0CLIM(i,j)=max(S0CLIM(i,j),0.)   ! pour etre au niveau des mers
      end do
   end do
   close(20)

   ! temperature : temperature a l'altitude de la calotte de ref (MAR)
   filin=trim(dirnameinp)//trim(filtr_t)
   write(6,*) filin
   open(20,file=filin)
   do j=1,ny
      do i=1,nx
         read(20,*) ti, tj, (Tm(i,j,mo),mo=1,12)
      end do
   end do
   close(20)

! test: on garde la temp parametree ->
!!$   do mo=1,12
!!$      Tm(:,:,mo)=Tm_surf(:,:,mo)
!!$   end do


   ! precipitation a l'altitude de la calotte de ref (MAR)
   filin=trim(dirnameinp)//trim(filtr_p)
   write(6,*) filin
   open(20,file=filin)
   do j=1,ny
      do i=1,nx
         read(20,*) ti, tj, (Pm(i,j,mo),mo=1,12)
      end do
   end do
   close(20)



   do j=1,ny
      do i=1,nx
         do mo=1,mois
            Tm_0(i,j,mo)=Tm(i,j,mo)      ! temp sur la calotte du GCM a l'instant initial
            
            zel=max(sealevel_2d(i,j),S0CLIM(i,j))
! on calcule une temperature au sol a partir de la temperature sur la calotte de ref du GCM


            if(zel.gt.3000) then
               Tsol(i,j,mo)=Tm(i,j,mo)-(lapserate(mo,4)*(zel-3000)+lapserate(mo,3)*1000+lapserate(mo,2)*1000+lapserate(mo,1)*1000)
            end if
            if ((zel.le.3000).and.(zel.gt.2000)) then
               Tsol(i,j,mo)=Tm(i,j,mo)-(lapserate(mo,3)*(zel-2000)+lapserate(mo,2)*1000+lapserate(mo,1)*1000)
            end if
            if ((zel.le.2000).and.(zel.gt.1000)) then
               Tsol(i,j,mo)=Tm(i,j,mo)-(lapserate(mo,2)*(zel-1000)+lapserate(mo,1)*1000)
            end if
            if (zel.le.1000) then
               Tsol(i,j,mo)=Tm(i,j,mo)-lapserate(mo,1)*zel
            end if

            !zel=max(sealevel_2d(i,j),S(i,j))    
            zel=max(0.,S(i,j))    
! puis on calcule la temperature sur la calotte S

            if(zel.gt.3000) then
               T_surf_ref=Tsol(i,j,mo)+(lapserate(mo,4)*(zel-3000)+lapserate(mo,3)*1000+lapserate(mo,2)*1000+lapserate(mo,1)*1000)
            end if
            if((zel.le.3000).and.(zel.gt.2000)) then
               T_surf_ref=Tsol(i,j,mo)+(lapserate(mo,3)*(zel-2000)+lapserate(mo,2)*1000+lapserate(mo,1)*1000)
            end if
            if ((zel.le.2000).and.(zel.gt.1000)) then
               T_surf_ref=Tsol(i,j,mo)+(lapserate(mo,2)*(zel-1000)+lapserate(mo,1)*1000)
            end if
            if (zel.le.1000) then
               T_surf_ref=Tsol(i,j,mo)+lapserate(mo,1)*zel
            end if

            if(retroac.eq.1) then      ! full retroaction acc.
               Pm_surf(i,j,mo)= Pm(i,j,mo)*exp(rappact*(T_surf_ref-Tm(i,j,mo)))
            else if (retroac.eq.0) then
               Pm_surf(i,j,mo)=Pm(i,j,mo)
            end if

            Tm_surf(i,j,mo)=T_surf_ref
            
         end do
      end do
   end do




Tm_surf(:,:,mois+1)=Tm_surf(:,:,1)
!Tm_0(:,:,mois+1)=Tm_0(:,:,1)
Pm_surf(:,:,mois+1)=Pm_surf(:,:,1)
Tsol(:,:,mois+1)=Tsol(:,:,1)

do j=1,ny
   do i=1,nx
      Tann(i,j)=(sum(Tm_surf(i,j,:))-Tm_surf(i,j,1))/12 ! janvier compte en double (modulo)
      Tjuly(i,j)=Tm_surf(i,j,7)
      !Tjuly(i,j)=(Tm_surf(i,j,6)+Tm_surf(i,j,7)+Tm_surf(i,j,8))/3 ! en fait un Tjja
   end do
end do

!!$do j=25,40
!!$do i=20,30
!!$print*, 'ecart annuel ', i, j, Tann(i,j),Tann(i,j)-(sum(Tm(i,j,:))-Tm(i,j,1))/12
!!$end do
!!$end do

print*, 'fini input_climat_ref'



end subroutine input_climat_ref




!------------------------------------------------------------------------------------------------------------
subroutine lect_lapserate_months  ! lapserates mensuels mais uniformes spatialement
! cette routine est susceptible de changer en fonction de la facon de lire les lapserates.
! pour le moment : on peut lire un fichier contenant les lapserates non homogènes spatiallement.


character (len=200) :: filtr_lr
integer :: k

real l1ann,l2ann,l3ann,l4ann
real l1july,l2july,l3july,l4july

!!$namelist/lapse_month/filtr_lr
!!$
!!$! lecture de la namelist
!!$! formats pour les ecritures dans 42
!!$428 format(A)
!!$
!!$rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
!!$read(num_param,lapse_month)
!!$
!!$write(num_rep_42,428)'!___________________________________________________________' 
!!$write(num_rep_42,428) '&lapse_month                                          ! module climat_perturb_mois_mod'
!!$write(num_rep_42,'(A,A)')   'filtr_t =', filtr_lr
!!$write(num_rep_42,*)'/'                      
!!$write(num_rep_42,428) '! lapse rates janvier -> decembre en deg/m'
!!$write(num_rep_42,*)
!!$
!!$
!!$! 3 colonnes de 12 valeurs (mois)
!!$! les 3 colonnes correspondent aux 3 tranches d'altitude
!!$! 0-1000m // 1000m-2000m // 2000m-++
!!$open(20,file=trim(dirnameinp)//trim(filtr_lr))
!!$do mo=1,mois
!!$      read(20,*) (lapserate(mo,k),k=1,3) ! en deg/m
!!$      ! test aurel pour supprimer la correction d'alti
!!$      !do k=1,3
!!$      !   lapserate(mo,k)=0.
!!$      !enddo
!!$      ! fin test aurel sur lapse
!!$end do
!!$close(20)

! le fichier n'ayant pas specialement des bonnes valeurs
! on utilise plutot les valeurs de Gherard 
! 4 tranches <1000, 1000-2000, 2000-3000, >3000 (l1, l2, l3, l4)

l1ann=-5.3/1000.   ! <1000
l2ann=-5.9/1000.   ! 1000-2000
l3ann=-8.1/1000.   ! 2000-3000
l4ann=-7.9/1000.   ! > 3000
! on suppose que le gradient de juillet est plus faible que l'annuel
! le 14% vient de Fausto et al. 2009
l1july=l1ann-l1ann*0.14
l2july=l2ann-l2ann*0.14
l3july=l3ann-l3ann*0.14
l4july=l4ann-l4ann*0.14


do mo=1,12
   lapserate(mo,1)=l1ann+(l1july-l1ann)*cos(2*pi*(mo+5)/12) 
   lapserate(mo,2)=l2ann+(l2july-l2ann)*cos(2*pi*(mo+5)/12) 
   lapserate(mo,3)=l3ann+(l3july-l3ann)*cos(2*pi*(mo+5)/12) 
   lapserate(mo,4)=l4ann+(l4july-l4ann)*cos(2*pi*(mo+5)/12)
end do


end subroutine lect_lapserate_months
!--------------------------------------------------------------------------------------------------------



end module lect_clim_act_greeneem_mois_lapsecouche