Index: /branches/iLoveclim/INPUT/ANTEIS1/.gmtcommands4
===================================================================
--- /branches/iLoveclim/INPUT/ANTEIS1/.gmtcommands4	(revision 51)
+++ /branches/iLoveclim/INPUT/ANTEIS1/.gmtcommands4	(revision 52)
@@ -1,2 +1,3 @@
 # GMT common arguments shelf
+-R1/141/1/141
 EOF
Index: /branches/iLoveclim/SOURCES/Ant40_files/bmelt-ant-regions_mod.f90
===================================================================
--- /branches/iLoveclim/SOURCES/Ant40_files/bmelt-ant-regions_mod.f90	(revision 51)
+++ /branches/iLoveclim/SOURCES/Ant40_files/bmelt-ant-regions_mod.f90	(revision 52)
@@ -14,298 +14,332 @@
 MODULE  BMELT_ANT_REGIONS         ! cat juillet 2005
 
-       USE module3D_phy
-
-
-implicit none
-
-REAL,dimension(nx,ny) ::  bmgrz       !< fusion basale a la grounding zone
-real,dimension(nx,ny) ::  bmshelf     !< fusion basale sous shelf
-
-REAL,dimension(nx,ny) :: dist_talu    !< distance du point au talu continental
-REAL,dimension(nx,ny) :: typeshelf    !< Type de shelf Ronne->1 Ross ->2 ....
-                                      !< utilises pour moduler la fusion sous le shelf
-integer, dimension(10) :: region      !< pour écrire dans le fichier param
-character(len=30),dimension(10) :: regname !< nom des régions 
-real :: bsupshelf
-integer :: grd                        !< pour une sortie 
-
-real :: bmelt_Ross,bmgrz_Ross         !< fusion basale Ross
-real :: bmelt_FRis,bmgrz_FRis         !< fusion basale filchner-Ronne ice shelf
-real :: bmelt_Amery,bmgrz_Amery       !< Amery ice shelves
-real :: bmelt_PIG,bmgrz_PIG           !< Pine Island glacier
-real :: bmelt_Pen,bmgrz_Pen           !< Petits ice shelves Peninsule
-real :: bmelt_other,bmgrz_other       !< Autre ice shelves
-real :: bmelt_talus,bmgrz_talus       !< Au dela du talus continental
+  USE module3D_phy
+
+
+  implicit none
+
+  REAL,dimension(nx,ny) ::  bmgrz       !< fusion basale a la grounding zone
+  real,dimension(nx,ny) ::  bmshelf     !< fusion basale sous shelf
+
+  REAL,dimension(nx,ny) :: dist_talu    !< distance du point au talu continental
+  REAL,dimension(nx,ny) :: typeshelf    !< Type de shelf Ronne->1 Ross ->2 ....
+  !< utilises pour moduler la fusion sous le shelf
+  integer, dimension(10) :: region      !< pour écrire dans le fichier param
+  character(len=30),dimension(10) :: regname !< nom des régions 
+  real :: bsupshelf
+  integer :: grd                        !< pour une sortie 
+
+  real :: bmelt_Ross,bmgrz_Ross         !< fusion basale Ross
+  real :: bmelt_FRis,bmgrz_FRis         !< fusion basale filchner-Ronne ice shelf
+  real :: bmelt_Amery,bmgrz_Amery       !< Amery ice shelves
+  real :: bmelt_PIG,bmgrz_PIG           !< Pine Island glacier
+  real :: bmelt_Pen,bmgrz_Pen           !< Petits ice shelves Peninsule
+  real :: bmelt_other,bmgrz_other       !< Autre ice shelves
+  real :: bmelt_talus,bmgrz_talus       !< Au dela du talus continental
+  real :: bmelt_coef                    !< Coef pour corriger les valeurs std de bmelt (lu dans namelist)
 
 CONTAINS
-!-------------------------------------------------------------------------------
-subroutine init_bmelt
-
-
-
-! Cette routine fait l'initialisation pour la fusion basale.
-! Elle est appelée par inputfile-vec-0.5.f90
-
-
-namelist/bmelt_ant_reg/bmelt_Ross,bmgrz_Ross,bmelt_FRis,bmgrz_FRis,  & 
-                       bmelt_Amery,bmgrz_Amery,bmelt_PIG,bmgrz_PIG,  &
-                       bmelt_Pen,bmgrz_Pen,bmelt_other,bmgrz_other,  &       
-                       bmelt_talus,bmgrz_talus      
-
-
-rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
-read(num_param,bmelt_ant_reg)
-
-
-!    ecriture dans le fichier parametres
-     write(num_rep_42,*)'fusion basale sous les ice shelves : bmelt-ant-regions_mod'
-     write(num_rep_42,*)'----------------------------------------------------------'
-
-! lecture du fichier contenant les distances au talu continental (m)
-      open(88,file=TRIM(DIRNAMEINP)//'distance_talu-40km.xy')
-
-      do j=1,ny
-         do i=1,nx
-            read(88,'(i3,1x,i3,1x,f10.2)') k,k,dist_talu(i,j)
-         enddo
-      enddo
-      close(88)
-
-! les bords de la grille sont mis a dist_talus=0 de force
-
-
-      dist_talu(1,:)=0.   !a gauche
-      dist_talu(2,:)=0.
-      
-      dist_talu(nx-1,:)=0. !  a droite 
-      dist_talu(nx,:)=0.  
-     
-      dist_talu(:,1)=0.    ! en bas   
-      dist_talu(:,2)=0.
-
-      dist_talu(:,ny-1)=0.    ! en haut
-      dist_talu(:,ny)=0.
-
-
-      debug_3D(:,:,32)=dist_talu(:,:)
-
-! lecture du fichier contenant les types de shelves
-!  Ronne-Flichner ->1, Ross -> 2 , Amery -> 3, 
-!  PIG-> 4, les petits shelves au dessus de PIG -> 5
-
-      open(88,file=TRIM(DIRNAMEINP)//'numer-ice-shelves-juil07.dat',status='OLD')
-      typeshelf(:,:)=100
-      do k=1,nx*ny
-        read(88,*,end=36) i,j,typeshelf(i,j)
-      end do
-36      close(88)
-
-
-
-  region(:)=0
-  regname(1)='Ronne-Fichner'
-  regname(2)='Ross'
-  regname(3)='Amery'
-  regname(4)='Pig'
-  regname(5)='Petits ice shelves peninsule'
-  regname(6)='autres ice shelves'
-  regname(7)='Au dela du talus continental'
- 
-bms_init:    do j=1,ny
-         do i=1,nx
-
-talus:    if (dist_talu(i,j).GT.5000.) then
-            if (nint(typeshelf(i,j)).eq.1) then ! Ronne-Filchner FRis
-                    bmshelf(i,j)=bmelt_FRis         
-                    bmgrz(i,j)=bmgrz_FRis         
-
-
-               if (region(1).eq.0) then          ! pour l'impression dans 42
-                  region(1)=1
-                  write(num_rep_42,80)regname(1),bmshelf(i,j),bmgrz(i,j)
-               endif
-
-80    format(a32,' bmshelf=',f8.4,'  bmgrz=',f8.4)
-
-            else  if (nint(typeshelf(i,j)).eq.2) then ! Ross
-               bmshelf(i,j)=bmelt_Ross  
-               bmgrz(i,j)=bmgrz_Ross      
-
-               if (region(2).eq.0) then
-                  region(2)=1
-                  write(num_rep_42,80)regname(2),bmshelf(i,j),bmgrz(i,j)
-               endif
-
-
-            else  if (nint(typeshelf(i,j)).eq.3) then ! Amery
-                    bmshelf(i,j)=bmelt_Amery   
-                    bmgrz(i,j)=bmgrz_Amery          
-               if (region(3).eq.0) then
-                  region(3)=1
-                  write(num_rep_42,80)regname(3),bmshelf(i,j),bmgrz(i,j)
-               endif
-
-
-            else  if (nint(typeshelf(i,j)).eq.4) then ! Pig
-                    bmshelf(i,j)=bmelt_PIG  
-                    bmgrz(i,j)=bmgrz_PIG       
-
-               if (region(4).eq.0) then
-                  region(4)=1
-                  write(num_rep_42,80)regname(4),bmshelf(i,j),bmgrz(i,j)
-               endif
-
-            else  if (nint(typeshelf(i,j)).eq.5) then ! petits shelves peninsule
-                    bmshelf(i,j)=bmelt_Pen     
-                    bmgrz(i,j)=bmgrz_Pen       
-
-               if (region(5).eq.0) then
-                  region(5)=1
-                  write(num_rep_42,80)regname(5),bmshelf(i,j),bmgrz(i,j)
-               endif
-
-
-            else !
-                    typeshelf(i,j)=6
-                    bmshelf(i,j)=bmelt_other     
-                    bmgrz(i,j)=bmgrz_other       
-
-               if (region(6).eq.0) then
-                  region(6)=1
-                  write(num_rep_42,80)regname(6),bmshelf(i,j),bmgrz(i,j)
-               endif
-
-
-            endif
-
-
-       else   ! au dela du talus continental
-
-                    bmshelf(i,j)=bmelt_talus  
-                    bmgrz(i,j)=bmgrz_talus   
- 
-               if (region(7).eq.0) then
-                  region(7)=1
-                  write(num_rep_42,80)regname(7),bmshelf(i,j),bmgrz(i,j)
-               endif
-
-
-            endif talus   ! fin du test sur dist_talus
-
-         enddo
-      enddo bms_init
-
-      debug_3D(:,:,33)=typeshelf(:,:)
-      debug_3D(:,:,34)=bmshelf(:,:)
-      return
-      end subroutine init_bmelt
-
-!________________________________________________________________________________
-subroutine bmeltshelf
-
-
-! cette routine calcule la fusion basale proprement dite
-
-integer :: ngr           ! nombre de voisins flottants
-real :: coef_talus       ! pour ne pas changer la fusion au dessus de l'ocean profond
+  !-------------------------------------------------------------------------------
+  subroutine init_bmelt
+
+
+    ! Cette routine fait l'initialisation pour la fusion basale.
+    ! Elle est appelée par inputfile-vec-0.5.f90
+
+
+    namelist/bmelt_ant_reg/bmelt_Ross,bmgrz_Ross,bmelt_FRis,bmgrz_FRis,  & 
+         bmelt_Amery,bmgrz_Amery,bmelt_PIG,bmgrz_PIG,  &
+         bmelt_Pen,bmgrz_Pen,bmelt_other,bmgrz_other,  &       
+         bmelt_talus,bmgrz_talus,bmelt_coef      
+
+
+    rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
+    read(num_param,bmelt_ant_reg)
+
+    ! formats pour les ecritures dans 42
+428 format(A)
+
+    write(num_rep_42,428)'!___________________________________________________________' 
+    write(num_rep_42,428) '&bmelt-ant-regions                                         '
+    write(num_rep_42,*)
+    write(num_rep_42,*) 'bmelt_Ross      = ',bmelt_Ross 
+    write(num_rep_42,*) 'bmgrz_Ross      = ',bmgrz_Ross
+    write(num_rep_42,*) 'bmelt_FRis      = ',bmelt_FRis 
+    write(num_rep_42,*) 'bmgrz_FRis      = ',bmgrz_FRis
+    write(num_rep_42,*) 'bmelt_Amery     = ',bmelt_Amery 
+    write(num_rep_42,*) 'bmgrz_Amery     = ',bmgrz_Amery
+    write(num_rep_42,*) 'bmelt_PIG       = ',bmelt_PIG 
+    write(num_rep_42,*) 'bmgrz_PIG       = ',bmgrz_PIG
+    write(num_rep_42,*) 'bmelt_Pen       = ',bmelt_Pen
+    write(num_rep_42,*) 'bmgrz_Pen       = ',bmgrz_Pen
+    write(num_rep_42,*) 'bmelt_other     = ',bmelt_other 
+    write(num_rep_42,*) 'bmgrz_other     = ',bmgrz_other
+    write(num_rep_42,*) 'bmelt_talus     = ',bmelt_talus 
+    write(num_rep_42,*) 'bmgrz_talus     = ',bmgrz_talus
+    write(num_rep_42,*) 'bmelt_coef      = ',bmelt_coef
+    write(num_rep_42,*)'/'                            
+    write(num_rep_42,428) '! bmelt_Ross & bmgrz_Ross    :  fonte basale secteur Ross'
+    write(num_rep_42,428) '! bmelt_FRis & bmgrz_FRis    :  fonte basale secteur Filchner Ronne'
+    write(num_rep_42,428) '! bmelt_Amery & bmgrz_Amery  :  fonte basale secteur Amery'
+    write(num_rep_42,428) '! bmelt_PIG & bmgrz_PIG      :  fonte basale secteur Pine Island'
+    write(num_rep_42,428) '! bmelt_Pen & bmgrz_Pen      :  fonte basale secteur Peninsule'
+    write(num_rep_42,428) '! bmelt_other & bmgrz_other  :  fonte basale autres secteurs'
+    write(num_rep_42,428) '! bmelt_talus & bmgrz_talus  :  fonte basale apres talus cont'
+    write(num_rep_42,428) '! bmelt_coef                 :  coef fonte (1 pour conserver val)'
+    write(num_rep_42,*)
+
+
+    !    ecriture dans le fichier parametres
+    write(num_rep_42,*)'fusion basale sous les ice shelves : bmelt-ant-regions_mod'
+    write(num_rep_42,*)'----------------------------------------------------------'
+
+
+    ! lecture du fichier contenant les distances au talu continental (m)
+    open(88,file=TRIM(DIRNAMEINP)//'distance_talu-40km.xy')
+
+    do j=1,ny
+       do i=1,nx
+          read(88,'(i3,1x,i3,1x,f10.2)') k,k,dist_talu(i,j)
+       enddo
+    enddo
+    close(88)
+
+    ! les bords de la grille sont mis a dist_talus=0 de force
+
+
+    dist_talu(1,:)=0.   !a gauche
+    dist_talu(2,:)=0.
+
+    dist_talu(nx-1,:)=0. !  a droite 
+    dist_talu(nx,:)=0.  
+
+    dist_talu(:,1)=0.    ! en bas   
+    dist_talu(:,2)=0.
+
+    dist_talu(:,ny-1)=0.    ! en haut
+    dist_talu(:,ny)=0.
+
+
+    debug_3D(:,:,32)=dist_talu(:,:)
+
+    ! lecture du fichier contenant les types de shelves Ronne-Flichner
+    !  ->1, Ross -> 2 , Amery -> 3, PIG-> 4, les petits shelves au dessus de PIG -> 5
+
+    open(88,file=TRIM(DIRNAMEINP)//'numer-ice-shelves-juil07.dat',status='OLD')
+    typeshelf(:,:)=100
+    do k=1,nx*ny
+       read(88,*,end=36) i,j,typeshelf(i,j)
+    end do
+36  close(88)
+
+
+
+    region(:)=0
+    regname(1)='Ronne-Fichner'
+    regname(2)='Ross'
+    regname(3)='Amery'
+    regname(4)='Pig'
+    regname(5)='Petits ice shelves peninsule'
+    regname(6)='autres ice shelves'
+    regname(7)='Au dela du talus continental'
+
+    bms_init:    do j=1,ny
+       do i=1,nx
+
+          talus:    if (dist_talu(i,j).GT.5000.) then
+             if (nint(typeshelf(i,j)).eq.1) then ! Ronne-Filchner FRis
+                bmshelf(i,j)=bmelt_FRis         
+                bmgrz(i,j)=bmgrz_FRis         
+
+
+                if (region(1).eq.0) then          ! pour l'impression dans 42
+                   region(1)=1
+                   write(num_rep_42,80)regname(1),bmshelf(i,j),bmgrz(i,j)
+                endif
+
+80              format(a32,' bmshelf=',f8.4,'  bmgrz=',f8.4)
+
+             else  if (nint(typeshelf(i,j)).eq.2) then ! Ross
+                bmshelf(i,j)=bmelt_Ross  
+                bmgrz(i,j)=bmgrz_Ross      
+
+                if (region(2).eq.0) then
+                   region(2)=1
+                   write(num_rep_42,80)regname(2),bmshelf(i,j),bmgrz(i,j)
+                endif
+
+
+             else  if (nint(typeshelf(i,j)).eq.3) then ! Amery
+                bmshelf(i,j)=bmelt_Amery   
+                bmgrz(i,j)=bmgrz_Amery          
+                if (region(3).eq.0) then
+                   region(3)=1
+                   write(num_rep_42,80)regname(3),bmshelf(i,j),bmgrz(i,j)
+                endif
+
+
+             else  if (nint(typeshelf(i,j)).eq.4) then ! Pig
+                bmshelf(i,j)=bmelt_PIG  
+                bmgrz(i,j)=bmgrz_PIG       
+
+                if (region(4).eq.0) then
+                   region(4)=1
+                   write(num_rep_42,80)regname(4),bmshelf(i,j),bmgrz(i,j)
+                endif
+
+             else  if (nint(typeshelf(i,j)).eq.5) then ! petits shelves peninsule
+                bmshelf(i,j)=bmelt_Pen     
+                bmgrz(i,j)=bmgrz_Pen       
+
+                if (region(5).eq.0) then
+                   region(5)=1
+                   write(num_rep_42,80)regname(5),bmshelf(i,j),bmgrz(i,j)
+                endif
+
+
+             else !
+                typeshelf(i,j)=6
+                bmshelf(i,j)=bmelt_other     
+                bmgrz(i,j)=bmgrz_other       
+
+                if (region(6).eq.0) then
+                   region(6)=1
+                   write(num_rep_42,80)regname(6),bmshelf(i,j),bmgrz(i,j)
+                endif
+
+
+             endif
+
+
+          else   ! au dela du talus continental
+
+             bmshelf(i,j)=bmelt_talus  
+             bmgrz(i,j)=bmgrz_talus   
+
+             if (region(7).eq.0) then
+                region(7)=1
+                write(num_rep_42,80)regname(7),bmshelf(i,j),bmgrz(i,j)
+             endif
+
+
+          endif talus   ! fin du test sur dist_talus
+
+       enddo
+    enddo bms_init
+
+    bmshelf(:,:)=bmshelf(:,:)*bmelt_coef
+    bmgrz(:,:)=bmgrz(:,:)*bmelt_coef
+    debug_3D(:,:,33)=typeshelf(:,:)
+    debug_3D(:,:,34)=bmshelf(:,:)
+    return
+  end subroutine init_bmelt
+
+  !________________________________________________________________________________
+  subroutine bmeltshelf
+
+
+    ! cette routine calcule la fusion basale proprement dite
+
+    integer :: ngr           ! nombre de voisins flottants
+    real :: coef_talus       ! pour ne pas changer la fusion au dessus de l'ocean profond
 
     do j=2,ny-1
-      do i=2,nx-1
- 
-talus_nochange : if (typeshelf(i,j).lt.10) then 
-           coef_talus = coefbmshelf
-        else
-           coef_talus = 1         ! pas de changement au dela du talus continental
-        endif talus_nochange
-
-
-shelf:    if (flot(i,j)) then    ! partie flottante
-
-            bmelt(i,j)=coef_talus*bmshelf(i,j)
-
-! bloc bsupshelf à enlever
-!                  if (time.gt.5000.) then
-!                     bmelt(i,j)=bmelt(i,j)+bsupshelf
-!                  endif    
-
-
-
-            if (fbm(i,j)) then 
-               bmelt(i,j)=coef_talus*bmgrz(i,j)
-            endif   
-
-!                  if (time.gt.5000.) then
-!                     bmelt(i,j)=bmelt(i,j)+bsupshelf
-!                  endif    
-
-! ATTENTION le bloc suivant sert a determiner la fusion basale d'equilibre
-! pour les shelves stationnaires
-
-  if (igrdline.eq.1) then
-     corrbmelt(i,j)=hdot(i,j)+bmelt(i,j)   ! le bmelt d'equilibre
-     debug_3D(i,j,28)=corrbmelt(i,j)
-  endif
-
-
-        else                   ! point posé, on compte le nombre de voisins flottants
-           ngr=0
-           if (flot(i+1,j)) ngr=ngr+1
-           if (flot(i-1,j)) ngr=ngr+1
-           if (flot(i,j+1)) ngr=ngr+1
-           if (flot(i,j-1)) ngr=ngr+1
-
-!   la fusion des points limites est une combinaison entre valeur posée et valeur flottante
-!   en fonction du nombre de points flottants
-
-           bmelt(i,j)= ngr/4.*bmgrz(i,j)*coef_talus+(1.-ngr/4.)*bmelt(i,j)
-         
-
-        endif shelf
-
-      end do 
+       do i=2,nx-1
+
+          talus_nochange : if (typeshelf(i,j).lt.10) then 
+             coef_talus = coefbmshelf
+          else
+             coef_talus = 1         ! pas de changement au dela du talus continental
+          endif talus_nochange
+
+
+          shelf:    if (flot(i,j)) then    ! partie flottante
+
+             bmelt(i,j)=coef_talus*bmshelf(i,j)
+
+             ! bloc bsupshelf à enlever
+             !                  if (time.gt.5000.) then
+             !                     bmelt(i,j)=bmelt(i,j)+bsupshelf
+             !                  endif    
+
+
+
+             if (fbm(i,j)) then 
+                bmelt(i,j)=coef_talus*bmgrz(i,j)
+             endif
+
+             !                  if (time.gt.5000.) then
+             !                     bmelt(i,j)=bmelt(i,j)+bsupshelf
+             !                  endif    
+
+             ! ATTENTION le bloc suivant sert a determiner la fusion basale d'equilibre
+             ! pour les shelves stationnaires
+
+             if (igrdline.eq.1) then
+                corrbmelt(i,j)=hdot(i,j)+bmelt(i,j)   ! le bmelt d'equilibre
+                debug_3D(i,j,28)=corrbmelt(i,j)
+             endif
+
+
+          else                   ! point posé, on compte le nombre de voisins flottants
+             ngr=0
+             if (flot(i+1,j)) ngr=ngr+1
+             if (flot(i-1,j)) ngr=ngr+1
+             if (flot(i,j+1)) ngr=ngr+1
+             if (flot(i,j-1)) ngr=ngr+1
+
+             !   la fusion des points limites est une combinaison entre valeur posée et valeur flottante
+             !   en fonction du nombre de points flottants
+
+             bmelt(i,j)= ngr/4.*bmgrz(i,j)*coef_talus+(1.-ngr/4.)*bmelt(i,j)
+
+
+          endif shelf
+
+       end do
     end do
 
-! les bords de la grille sont mis a dist_talus=0 de force
-! attention, pourrait poser des problemes avec AGRIF
-
-      bmelt(1,:)=bmelt_talus   !a gauche
-      bmelt(2,:)=bmelt_talus
-      
-      bmelt(nx-1,:)=bmelt_talus !  a droite 
-      bmelt(nx,:)=bmelt_talus  
-     
-      bmelt(:,1)=bmelt_talus    ! en bas   
-      bmelt(:,2)=bmelt_talus
-
-      bmelt(:,ny-1)=bmelt_talus    ! en haut
-      bmelt(:,ny)=bmelt_talus
-
-
-
-if (igrdline.eq.1) then
-    bmelt(:,:)=0.      ! hdot donne alors le -bmelt
-endif
-
-! ecriture des bmelt pour en faire des statistiques
-!--------------------------------------------------
-!open(144,file='bmelt-test')
-!do j=1,ny
-!   do i=1,nx
-!      if (fbm(i,j)) then
-!         grd=1
-!      else
-!         grd=0
-!      endif
-!      if ((flot(i,j)).and.(H(i,j).gt.1.)) then
-!         write(144,999) i,j,grd,nint(typeshelf(i,j)), & 
-!              bmelt(i,j),bmgrz(i,j),bmshelf(i,j),hdot(i,j),dist_talu(i,j)
-!      endif
-!   end do
-!end do
-!999 format(4(i0,2x),5(e12.3,2x))
-!close(144)
+    ! les bords de la grille sont mis a dist_talus=0 de force
+    ! attention, pourrait poser des problemes avec AGRIF
+
+    bmelt(1,:)=bmelt_talus   !a gauche
+    bmelt(2,:)=bmelt_talus
+
+    bmelt(nx-1,:)=bmelt_talus !  a droite 
+    bmelt(nx,:)=bmelt_talus  
+
+    bmelt(:,1)=bmelt_talus    ! en bas   
+    bmelt(:,2)=bmelt_talus
+
+    bmelt(:,ny-1)=bmelt_talus    ! en haut
+    bmelt(:,ny)=bmelt_talus
+
+
+
+    if (igrdline.eq.1) then
+       bmelt(:,:)=0.      ! hdot donne alors le -bmelt
+    endif
+
+    ! ecriture des bmelt pour en faire des statistiques
+    !--------------------------------------------------
+    !open(144,file='bmelt-test')
+    !do j=1,ny
+    !   do i=1,nx
+    !      if (fbm(i,j)) then
+    !         grd=1
+    !      else
+    !         grd=0
+    !      endif
+    !      if ((flot(i,j)).and.(H(i,j).gt.1.)) then
+    !         write(144,999) i,j,grd,nint(typeshelf(i,j)), & 
+    !              bmelt(i,j),bmgrz(i,j),bmshelf(i,j),hdot(i,j),dist_talu(i,j)
+    !      endif
+    !   end do
+    !end do
+    !999 format(4(i0,2x),5(e12.3,2x))
+    !close(144)
 
     return
-    end subroutine  bmeltshelf
+  end subroutine  bmeltshelf
 
 
Index: /branches/iLoveclim/SOURCES/Ant40_files/lect-anteis_mod.f90
===================================================================
--- /branches/iLoveclim/SOURCES/Ant40_files/lect-anteis_mod.f90	(revision 51)
+++ /branches/iLoveclim/SOURCES/Ant40_files/lect-anteis_mod.f90	(revision 52)
@@ -242,7 +242,5 @@
 !    close(88)
 
-    call lect_input(1,'ghf',1,ghf,ghf_fich,file_ncdf)    ! socle
-
-
+    call lect_input(1,'ghf',1,ghf,ghf_fich,file_ncdf)
 
     ! pour passer les flux des mW/m2 au J/m2/an      
Index: /branches/iLoveclim/SOURCES/Fichiers-parametres/anteis1_LISTE-VAR-NETCDF.dat
===================================================================
--- /branches/iLoveclim/SOURCES/Fichiers-parametres/anteis1_LISTE-VAR-NETCDF.dat	(revision 51)
+++ /branches/iLoveclim/SOURCES/Fichiers-parametres/anteis1_LISTE-VAR-NETCDF.dat	(revision 52)
@@ -118,5 +118,5 @@
 Acc		
 11	
-1	4	1				
+1	1	1				
 -----------------------------------------------------------
 Abl		
Index: /branches/iLoveclim/SOURCES/GrIce2sea_files/climat_GrIce2sea_years_perturb_mod.f90
===================================================================
--- /branches/iLoveclim/SOURCES/GrIce2sea_files/climat_GrIce2sea_years_perturb_mod.f90	(revision 51)
+++ /branches/iLoveclim/SOURCES/GrIce2sea_files/climat_GrIce2sea_years_perturb_mod.f90	(revision 52)
@@ -100,5 +100,5 @@
 
   namelist/clim_smb_T_gen/smb_file,coef_smb_unit,temp_annual_file
-
+  namelist/clim_snap/nb_snap,time_depart_snaps,ecart_snap,file_smb_snap,massb_time
 
 428 format(A)
@@ -118,5 +118,5 @@
   smb_file  = trim(dirnameinp)//trim(smb_file)
 
-  call Read_Ncdf_var('z',smb_file,tab)
+  call Read_Ncdf_var('smb',smb_file,tab)
 
 !  call lect_input(3,'smb',1,bm,smb_file,trim(dirnameinp)//trim(runname)//'.nc')
@@ -145,5 +145,5 @@
   temp_annual_file = trim(dirnameinp)//trim(temp_annual_file)
 
-  call Read_Ncdf_var('z',temp_annual_file,tab)
+  call Read_Ncdf_var('Tann',temp_annual_file,tab)
   Tann(:,:)  = tab(:,:)
 !  call lect_input(3,'Tann',1,Tann,temp_annual_file,trim(dirnameinp)//trim(runname)//'.nc')
@@ -155,10 +155,4 @@
   ta0(:,:)   = Tann(:,:)
   Tjuly(:,:) = Tann(:,:)
-
-
-
-!lecture namelist 
-
-  namelist/clim_snap/nb_snap,time_depart_snaps,ecart_snap,file_smb_snap,massb_time
 
   rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
Index: /branches/iLoveclim/SOURCES/Hemin40_files/lect-hemin40_mod.f90
===================================================================
--- /branches/iLoveclim/SOURCES/Hemin40_files/lect-hemin40_mod.f90	(revision 51)
+++ /branches/iLoveclim/SOURCES/Hemin40_files/lect-hemin40_mod.f90	(revision 52)
@@ -2,9 +2,16 @@
 
   use module3D_phy
+  use interface_input
+  use io_netcdf_grisli
+
+  implicit none
   
-    character(len=50) :: FILE1
-    character(len=50) :: FILE2
-    character(len=80) :: filin
-    real, dimension(nx,ny,5) :: bidon          ! pour l'appel a courbure
+  character(len=100) :: topo_dep       ! Topo de départ
+  character(len=100) :: topo_ref       ! Topo de référence
+  character(len=100) :: grid_topo      ! fichier grille
+  character(len=100) :: ghf_fich       ! fichier grille
+  character(len=80) :: filin
+  real, dimension(nx,ny,5) :: bidon          ! pour l'appel a courbure
+  character(len=100) :: file_ncdf      !< fichier netcdf issue des fichiers .dat
 
 contains
@@ -12,18 +19,36 @@
 subroutine input_topo
 
-namelist/topo_file/file1,file2
-rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
-read(num_param,topo_file)
+  integer :: ios
+
+  namelist/topo_file/topo_ref,topo_dep,grid_topo,ghf_fich
+  rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
+  read(num_param,topo_file)
 ! formats pour les ecritures dans 42
 428 format(A)
-write(num_rep_42,428)'!___________________________________________________________' 
-write(num_rep_42,428) '&topo_file                                  ! nom du bloc '
-write(num_rep_42,*)
-write(num_rep_42,*) 'file1 = ', file1
-write(num_rep_42,*) 'file2 = ', file2
-write(num_rep_42,*)'/' 
-write(num_rep_42,428) '! file1 : topo de depart'
-write(num_rep_42,428) '! file2 : topo de reference'             
-write(num_rep_42,*)
+  write(num_rep_42,428)'!___________________________________________________________' 
+  write(num_rep_42,428) '&topo_file                                  ! nom du bloc '
+  write(num_rep_42,*)
+  write(num_rep_42,'(A,A)') 'topo_ref = ', topo_ref
+  write(num_rep_42,'(A,A)') 'topo_dep = ', topo_dep
+  write(num_rep_42,'(A,A)') 'grid_topo =', grid_topo
+  write(num_rep_42,'(A,A)') 'ghf_fich = ', ghf_fich
+  write(num_rep_42,*)'/'                      
+  write(num_rep_42,428) '! topo_ref= topo ref isostasie'
+  write(num_rep_42,428) '! topo_dep= topo de depart'
+  write(num_rep_42,428) '! grid_topo : fichier i,j,x,y,lon,lat'
+  write(num_rep_42,428) '! ghf_fich  : fichier flux geothermique'  
+  write(num_rep_42,*)
+
+  topo_ref=trim(dirnameinp)//trim(topo_ref) 
+  topo_dep=trim(dirnameinp)//trim(topo_dep) 
+  grid_topo=trim(dirnameinp)//trim(grid_topo)
+  ghf_fich=trim(dirnameinp)//trim(ghf_fich)
+
+!write(num_rep_42,*) 'file1 = ', file1
+!write(num_rep_42,*) 'file2 = ', file2
+!write(num_rep_42,*)'/' 
+!write(num_rep_42,428) '! file1 : topo de depart'
+!write(num_rep_42,428) '! file2 : topo de reference'             
+!write(num_rep_42,*)
 
 !====================================== La reponse est 42 ===========
@@ -41,40 +66,51 @@
 
       
-! lecture adaptee aux fichiers intercomparaison EISMINT
-       nxx=nx
-       nyy=ny
+!!$! lecture adaptee aux fichiers intercomparaison EISMINT
+!!$       nxx=nx
+!!$       nyy=ny
 
-!     lecture de la topo actuelle
-!     ---------------------------
-     open (20,file=TRIM(DIRNAMEINP)//file2,status='old')
-       
-     read(20,'(A80)') TITRE
-     read(20,*) NI,NJ,NXX,NYY,STEP
-     read(20,*)
-         do J=1,ny 
-          do I=1,nx
-             read (20,*)  S0(I,J),H0(I,J),BSOC0(I,J)
-             S0(i,j)=max(S0(i,j),0.)
-          end do
-        end do
-     close(20)
+!!$!     lecture de la topo actuelle
+!!$!     ---------------------------
+!!$     open (20,file=TRIM(DIRNAMEINP)//file2,status='old')
+!!$       
+!!$     read(20,'(A80)') TITRE
+!!$     read(20,*) NI,NJ,NXX,NYY,STEP
+!!$     read(20,*)
+!!$         do J=1,ny 
+!!$          do I=1,nx
+!!$             read (20,*)  S0(I,J),H0(I,J),BSOC0(I,J)
+!!$             S0(i,j)=max(S0(i,j),0.)
+!!$          end do
+!!$        end do
+!!$     close(20)
+!!$
+!!$     
+!!$!     lecture de la topo de depart
+!!$!     ---------------------------
+!!$     open (20,file=TRIM(DIRNAMEINP)//file1,status='old')
+!!$!         open (20,file='../INPUT-DATA/hemin.g50')
+!!$     read(20,'(A80)') TITRE
+!!$     read(20,*) NI,NJ,NXX,NYY,STEP
+!!$     read(20,*)
+!!$         do J=1,ny 
+!!$          do I=1,nx
+!!$             read (20,*) S(I,J),H(I,J),BSOC(I,J)
+!!$          end do
+!!$        end do
+!!$     close(20)
 
-     
-!     lecture de la topo de depart
-!     ---------------------------
-     open (20,file=TRIM(DIRNAMEINP)//file1,status='old')
-!         open (20,file='../INPUT-DATA/hemin.g50')
-     read(20,'(A80)') TITRE
-     read(20,*) NI,NJ,NXX,NYY,STEP
-     read(20,*)
-         do J=1,ny 
-          do I=1,nx
-             read (20,*) S(I,J),H(I,J),BSOC(I,J)
-          end do
-        end do
-     close(20)
+
+!     lecture de la topo de référence
+    call lect_input(1,'Bsoc',1,Bsoc0,topo_ref,file_ncdf)    ! socle
+    call lect_input(1,'S',1,S0,topo_ref,file_ncdf)          ! surface
+    call lect_input(1,'H',1,H0,topo_ref,file_ncdf)          ! epaisseur
+
+!     lecture de la topo de départ
+    call lect_input(1,'Bsoc',1,Bsoc,topo_dep,file_ncdf)    ! socle
+    call lect_input(1,'S',1,S,topo_dep,file_ncdf)          ! surface
+    call lect_input(1,'H',1,H,topo_dep,file_ncdf)          ! epaisseur
+
 
 ! calcul des courbures du socle
-
      call courbure(nx,ny,dx,Bsoc,bidon(:,:,1),bidon(:,:,2),bidon(:,:,3), &
           bidon(:,:,4),socle_cry,bidon(:,:,5))
@@ -83,10 +119,10 @@
 ! lecture des coordonnées geographiques
 
-    filin=TRIM(DIRNAMEINP)//'coord-nord-40km.dat'
+!    filin=TRIM(DIRNAMEINP)//grid_topo
 
 ! les coordonnees sont calculees en °dec avec GMT,
 ! les longitudes sont comprises entre -180 et +180 (negative a l'Ouest de
 ! Greenwich et positive a l'Est)
-    open(unit=2004,file=filin,iostat=ios)
+    open(unit=2004,file=grid_topo,iostat=ios)
     do k=1,nx*ny
        read(2004,*) i,j,XCC(i,j),YCC(i,j),XLONG(i,j),YLAT(i,j)
@@ -101,13 +137,15 @@
                          
 ! lecture du flux geothermique de Shapiro
-    open(88,file=TRIM(DIRNAMEINP)//'ijphi_hemin40.dat')
+!    open(88,file=TRIM(DIRNAMEINP)//'ijphi_hemin40.dat')
+!    open(88,file=ghf_fich)
 
-    write(42,*) 'flux geothermique Shapiro : ',TRIM(DIRNAMEINP)//'ijphi_hemin40.dat'
+!!$    write(42,*) 'flux geothermique Shapiro : ',TRIM(DIRNAMEINP)//'ijphi_hemin40.dat'
+!!$    do k=1,nx*ny
+!!$       read(88,*) i,j,ghf(i,j)
+!!$    end do
+!!$    close(88)
 
-    do k=1,nx*ny
-       read(88,*) i,j,ghf(i,j)
-!        print*, i,j,ghf(i,j)
-    end do
-    close(88)
+    call lect_input(1,'ghf',1,ghf,ghf_fich,file_ncdf)
+
 ! pour passer les flux des mW/m2 au J/m2/an      
     ghf(:,:)=-SECYEAR/1000.*ghf(:,:)
Index: /branches/iLoveclim/SOURCES/Makefile.grisli-gfortran.inc
===================================================================
--- /branches/iLoveclim/SOURCES/Makefile.grisli-gfortran.inc	(revision 51)
+++ /branches/iLoveclim/SOURCES/Makefile.grisli-gfortran.inc	(revision 52)
@@ -92,5 +92,4 @@
         dragging_hwatermax_0.2_mod.o dragging_calc_beta_mod.o  \
         dragging_hwat-contigu_mod.o dragging_hwat_contmaj_mod.o \
-        prop-thermiques_mod.o  \
         calving_frange.o no_calving.o no_lakes.o \
         out_profile_mod.o printtable_mod.o mix-SIA-L1_mod.o \
@@ -1616,5 +1615,5 @@
 	$(toy_recul)  \
 	$(mod_communs_RGL) \
-	$(mod_clim_perturb)  \
+	$(mod_clim_tof)  \
 	$(mod_tracers) \
 	$(mod_ell) $(Liste_GrIce2sea) \
@@ -1628,5 +1627,5 @@
 	$(toy_recul)  \
 	$(mod_communs_RGL) \
-	$(mod_clim_perturb)  \
+	$(mod_clim_tof)  \
 	$(mod_tracers) \
 	$(mod_ell) $(Liste_GrIce2sea) \
Index: /branches/iLoveclim/SOURCES/Makefile.grisli.inc
===================================================================
--- /branches/iLoveclim/SOURCES/Makefile.grisli.inc	(revision 51)
+++ /branches/iLoveclim/SOURCES/Makefile.grisli.inc	(revision 52)
@@ -91,5 +91,4 @@
         dragging_hwatermax_0.2_mod.o dragging_calc_beta_mod.o  \
         dragging_hwat-contigu_mod.o dragging_hwat_contmaj_mod.o \
-        prop-thermiques_mod.o  \
         calving_frange.o no_calving.o no_lakes.o \
         out_profile_mod.o printtable_mod.o mix-SIA-L1_mod.o \
@@ -1617,5 +1616,5 @@
 	$(toy_recul)  \
 	$(mod_communs_RGL) \
-	$(mod_clim_perturb)  \
+	$(mod_clim_tof)  \
 	$(mod_tracers) \
 	$(mod_ell) $(Liste_GrIce2sea) \
@@ -1629,5 +1628,5 @@
 	$(toy_recul)  \
 	$(mod_communs_RGL) \
-	$(mod_clim_perturb)  \
+	$(mod_clim_tof)  \
 	$(mod_tracers) \
 	$(mod_ell) $(Liste_GrIce2sea) \
Index: /branches/iLoveclim/SOURCES/Netcdf-routines/io_netcdf_GRISLI.f90
===================================================================
--- /branches/iLoveclim/SOURCES/Netcdf-routines/io_netcdf_GRISLI.f90	(revision 51)
+++ /branches/iLoveclim/SOURCES/Netcdf-routines/io_netcdf_GRISLI.f90	(revision 52)
@@ -191,4 +191,8 @@
     !     
     status = nf90_inq_varid(ncid,varname,varid)
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif
     status=nf90_inquire_variable(ncid,varid,dimids=dimIDS)
     status=nf90_inquire_dimension(ncid,dimIDS(1),len=dim1)
@@ -249,4 +253,8 @@
     !      
     status = nf90_inq_varid(ncid,varname,varid)
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif
     status=nf90_inquire_variable(ncid,varid,dimids=dimIDS)
     !                
@@ -294,4 +302,8 @@
     !     
     status = nf90_inq_varid(ncid,varname,varid)
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif
     status=nf90_inquire_variable(ncid,varid,dimids=dimIDS)
     status=nf90_inquire_dimension(ncid,dimIDS(1),len=dim1)
@@ -346,4 +358,8 @@
     !     
     status = nf90_inq_varid(ncid,varname,varid)
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif
     status=nf90_inquire_variable(ncid,varid,dimids=dimIDS)
     status=nf90_inquire_dimension(ncid,dimIDS(1),len=dim1)
@@ -394,4 +410,8 @@
     !     
     status = nf90_inq_varid(ncid,varname,varid)
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif
     status=nf90_inquire_variable(ncid,varid,dimids=dimID)
     status=nf90_inquire_dimension(ncid,dimID(1),len=dim1)          
@@ -438,4 +458,8 @@
     !     
     status = nf90_inq_varid(ncid,varname,varid)
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif
     status=nf90_inquire_variable(ncid,varid,dimids=dimIDS)
     status=nf90_inquire_dimension(ncid,dimIDS(1),len=dim1)
@@ -484,4 +508,8 @@
     !     
     status = nf90_inq_varid(ncid,varname,varid)
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif
     status=nf90_inquire_variable(ncid,varid,dimids=dimIDS)
     status=nf90_inquire_dimension(ncid,dimIDS(1),len=dim1)
@@ -532,4 +560,8 @@
     !     
     status = nf90_inq_varid(ncid,varname,varid)
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif
     status=nf90_inquire_variable(ncid,varid,dimids=dimIDS)
     status=nf90_inquire_dimension(ncid,dimIDS(1),len=dim1)
@@ -584,5 +616,8 @@
     !     
     status = nf90_inq_varid(ncid,varname,varid)
-    !        
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif        
     status=nf90_inquire_variable(ncid,varid,dimids=dimIDS)
     status=nf90_inquire_dimension(ncid,dimIDS(1),len=dim1)
@@ -638,5 +673,8 @@
     !     
     status = nf90_inq_varid(ncid,varname,varid)
-    !        
+    if (status/=nf90_noerr) then
+       write(*,*)"pb de nom de variable :", varname, file
+       stop
+    endif        
     status=nf90_inquire_variable(ncid,varid,dimids=dimIDS)
     status=nf90_inquire_dimension(ncid,dimIDS(1),len=dim1)
Index: /branches/iLoveclim/SOURCES/steps_time_loop.f90
===================================================================
--- /branches/iLoveclim/SOURCES/steps_time_loop.f90	(revision 51)
+++ /branches/iLoveclim/SOURCES/steps_time_loop.f90	(revision 52)
@@ -15,9 +15,7 @@
 
   use module3d_phy
-
   use module_choix ! module de choix du type de run
-  !  module_choix donne acces a tous les modules
-  !  de declaration des packages
-
+                   !  module_choix donne acces a tous les modules
+                   !  de declaration des packages
   use sorties_ncdf_grisli
   use flottab_mod
Index: /branches/iLoveclim/SOURCES/steps_time_loop_avec_iterbeta.f90
===================================================================
--- /branches/iLoveclim/SOURCES/steps_time_loop_avec_iterbeta.f90	(revision 51)
+++ /branches/iLoveclim/SOURCES/steps_time_loop_avec_iterbeta.f90	(revision 52)
@@ -15,21 +15,11 @@
 
   use module3d_phy
-
-  use geom_type
-  use temperature_type
-  use ice_flow_type
-  use mask_flgz_type
-  use deformation_type
-  use autre_pr_temp_type
-
   use module_choix ! module de choix du type de run
                    !  module_choix donne acces a tous les modules
                    !  de declaration des packages
-
+  use icetempmod
   use sorties_ncdf_grisli
   use flottab_mod
-  use interface_icetempmod
   use diagno_mod 
-
 
 
@@ -205,5 +195,5 @@
   !  module_choix donne acces a tous les modules
   !  de declaration des packages
-  use interface_icetempmod
+  use icetempmod
   use sorties_ncdf_grisli
   use flottab_mod
@@ -211,12 +201,4 @@
   use resolmeca_SIA_L1
 !  use track_debug 
-
-  use geom_type
-  use temperature_type
-  use ice_flow_type
-  use mask_flgz_type
-  use deformation_type
-  use autre_pr_temp_type
-  use step_temp_type
 
 
@@ -309,29 +291,6 @@
 
      calc_temp:     if ((nt.gt.2).and.(geoplace(1:5).ne.'mism3')) then
-        ! Hassine ca va etre modifier bien sur apres: l appel sera beaucoup moin long que ca 
-        call Update_geom(Geom_g)               
-        call Update_Temperature(Temperature_g)
-        call Update_Ice_flow(Ice_flow_g)
-        call Update_Mask_flgz(Mask_flgz_g)
-        call Update_Deformation(Deformation_g)
-        call Update_autre_pr_temp(Autre_pr_temp_g)
-        call Update_step_temp(step_temp_g)
-
         if (itracebug.eq.1)  call tracebug('avant appel icetemp')
-
-        call icetemp(geom_g,temperature_g, &
-             Ice_flow_g,Mask_flgz_g,&
-             deformation_g,Autre_pr_temp_g,&
-             step_temp_g,num_tracebug,debug_3D)
-
-        ! update arrays T, Tpmp
-        T= Temperature_g%Temperature 
-        Tpmp= Temperature_g%Temp_melting
-        Tbdot= Autre_pr_temp_g%Tbdot
-        Bmelt= Autre_pr_temp_g%BMELT
-        Ibase=Autre_pr_temp_g%IBASE
-        Phid=Autre_pr_temp_g%Phid
-
-        !      write(6,*) time,'call icetemp'
+        call icetemp
         ! subroutines pour le calcul de la fusion basale
         call bmeltshelf
