Index: trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/ini_archive.F
===================================================================
--- trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/ini_archive.F	(revision 1644)
+++ trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/ini_archive.F	(revision 1647)
@@ -36,5 +36,4 @@
 
       USE comsoil_h
-      USE slab_ice_h, only: noceanmx
 !      use control_mod
       USE comvert_mod, ONLY: ap,bp,aps,bps,presnivs,pseudoalt
@@ -81,5 +80,4 @@
       INTEGER idim_tim
       INTEGER idim_nsoilmx ! "subsurface_layers" dimension ID #
-      INTEGER idim_noceanmx ! "ocean_layers" dimension ID #
       INTEGER nid,nvarid
       real sig_s(llm),s(llm)
@@ -164,5 +162,4 @@
       ierr = NF_DEF_DIM (nid, "altitude", llm, idim_llm)
       ierr = NF_DEF_DIM (nid,"subsurface_layers",nsoilmx,idim_nsoilmx)
-      ierr = NF_DEF_DIM (nid,"ocean_layers",noceanmx,idim_noceanmx)
 
       ierr = NF_DEF_DIM (nid,"index", length, idim_index)
Index: trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/iniphysiq_mod.F90
===================================================================
--- trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/iniphysiq_mod.F90	(revision 1644)
+++ trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/iniphysiq_mod.F90	(revision 1647)
@@ -11,5 +11,4 @@
 
 use control_mod, only: nday 
-use surf_heat_transp_mod, only: ini_surf_heat_transp
 use infotrac, only : nqtot ! number of advected tracers
 use planete_mod, only: ini_planete_mod
@@ -53,6 +52,4 @@
 integer,intent(in) :: iflag_phys ! type of physics to be called
 
-logical :: ok_slab_ocean
-
   ! the common part for all planetary physics
   !------------------------------------------
@@ -79,13 +76,4 @@
 call ini_planete_mod(nlayer,preff,ap,bp)
 
-! for slab ocean, copy over some arrays
-ok_slab_ocean=.false. ! default value
-call getin_p("ok_slab_ocean",ok_slab_ocean)
-if (ok_slab_ocean) then
-  call ini_surf_heat_transp(ip1jm,ip1jmp1,unsairez,fext,unsaire, &
-                            cu,cuvsurcv,cv,cvusurcu,aire,apoln,apols, &
-                            aireu,airev)
-endif
-
 call inifis(klon_omp,nlayer,nqtot,pdayref,punjours,nday,ptimestep, &
             latitude,longitude,cell_area,prad,pg,pr,pcpp)
Index: trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/lect_start_archive.F
===================================================================
--- trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/lect_start_archive.F	(revision 1644)
+++ trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/lect_start_archive.F	(revision 1647)
@@ -2,15 +2,11 @@
      &     date,tsurf,tsoil,emis,q2,
      &     t,ucov,vcov,ps,h,phisold_newgrid,
-     &     q,qsurf,surfith,nid,
-     &     rnat,pctsrf_sic,tslab,tsea_ice,sea_ice)
+     &     q,qsurf,surfith,nid)
 
 !      USE surfdat_h
       USE comsoil_h, ONLY: nsoilmx, layer, mlayer, volcapa, inertiedat
-      USE tracer_h, ONLY: igcm_co2_ice
       USE infotrac, ONLY: tname, nqtot
-      USE slab_ice_h, ONLY: noceanmx
 !      USE control_mod
 ! to use  'getin'
-      USE callkeys_mod, ONLY: ok_slab_ocean
       USE comvert_mod, ONLY: ap,bp,aps,bps,preff
       USE comconst_mod, ONLY: kappa,g,pi
@@ -101,10 +97,6 @@
       REAL tsurf(ngrid) ! surface temperature
       REAL tsoil(ngrid,nsoilmx) ! soil temperature
-      REAL co2ice(ngrid) ! CO2 ice layer
       REAL emis(ngrid)
       REAL q2(ngrid,llm+1),qsurf(ngrid,nqtot)
-      REAL tslab(ngrid,noceanmx)
-      REAL rnat(ngrid),pctsrf_sic(ngrid)
-      REAL tsea_ice(ngrid),sea_ice(ngrid)
 c     REAL phisfi(ngrid)
 
@@ -126,12 +118,8 @@
       real tsurfS(iip1,jjp1),tsoilS(iip1,jjp1,nsoilmx)
       real inertiedatS(iip1,jjp1,nsoilmx)
-      real co2iceS(iip1,jjp1)
       real emisS(iip1,jjp1)
       REAL q2S(iip1,jjp1,llm+1),qsurfS(iip1,jjp1,nqtot)
-      real tslabS(iip1,jjp1,noceanmx)
-      real pctsrf_sicS(iip1,jjp1),tsea_iceS(iip1,jjp1)
-      real rnatS(iip1,jjp1), sea_iceS(iip1,jjp1)
-
-      real ptotal, co2icetotal
+      
+      real ptotal
 
 c Var intermediaires : vent naturel, mais pas coord scalaire
@@ -156,16 +144,11 @@
       real, dimension(:,:,:), allocatable :: inertiedatoldnew
       real, dimension(:,:), allocatable :: psold,phisold
-      real, dimension(:,:), allocatable :: co2iceold
       real, dimension(:,:), allocatable :: tsurfold
       real, dimension(:,:), allocatable :: emisold
       real, dimension(:,:,:,:), allocatable :: qold
-      real, dimension(:,:,:), allocatable :: tslabold
-      real, dimension(:,:), allocatable :: rnatold,pctsrf_sicold
-      real, dimension(:,:), allocatable :: tsea_iceold,sea_iceold
-
 
       real tab_cntrl(100)
 
-      real ptotalold, co2icetotalold
+      real ptotalold
 
       logical :: olddepthdef=.false. ! flag
@@ -300,5 +283,4 @@
       allocate(phisold(imold+1,jmold+1))
       allocate(qold(imold+1,jmold+1,lmold,nqtot))
-      allocate(co2iceold(imold+1,jmold+1))
       allocate(tsurfold(imold+1,jmold+1))
       allocate(emisold(imold+1,jmold+1))
@@ -313,9 +295,4 @@
       allocate(mlayerold(nsoilold))
       allocate(qsurfold(imold+1,jmold+1,nqtot))
-      allocate(tslabold(imold+1,jmold+1,noceanmx))
-      allocate(rnatold(imold+1,jmold+1))
-      allocate(pctsrf_sicold(imold+1,jmold+1))
-      allocate(tsea_iceold(imold+1,jmold+1))
-      allocate(sea_iceold(imold+1,jmold+1))
 
       allocate(var (imold+1,jmold+1,llm))
@@ -332,5 +309,5 @@
 C-----------------------------------------------------------------------
 c 3.1. Lecture du tableau des parametres du run 
-c     (pour  la lecture ulterieure de "ptotalold" et "co2icetotalold")
+c     (pour  la lecture ulterieure de "ptotalold")
 c-----------------------------------------------------------------------
 c
@@ -556,8 +533,7 @@
 
 C-----------------------------------------------------------------------
-c   lecture de "ptotalold" et "co2icetotalold"
+c   lecture de "ptotalold"
 c-----------------------------------------------------------------------
       ptotalold = tab_cntrl(tab0+49)
-      co2icetotalold = tab_cntrl(tab0+50)
  
 c-----------------------------------------------------------------------
@@ -636,5 +612,5 @@
 
 c-----------------------------------------------------------------------
-c 5.1 Lecture des champs 2D (co2ice, emis,ps,tsurf,Tg[10], qsurf)
+c 5.1 Lecture des champs 2D (emis,ps,tsurf,Tg[10], qsurf)
 c-----------------------------------------------------------------------
  
@@ -642,21 +618,5 @@
       count=(/imold+1,jmold+1,1,0/)
        
-! CO2ice is now in qsurf(igcm_co2_ice) ...
-!      ierr = NF_INQ_VARID (nid, "co2ice", nvarid)
-!      IF (ierr .NE. NF_NOERR) THEN
-!         PRINT*, "lect_start_archive: Le champ <co2ice> est absent"
-!         CALL abort
-!      ENDIF
-!#ifdef NC_DOUBLE
-!      ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,count,co2iceold)
-!#else
-!      ierr = NF_GET_VARA_REAL(nid, nvarid,start,count,co2iceold)
-!#endif
-!      IF (ierr .NE. NF_NOERR) THEN
-!         PRINT*, "lect_start_archive: Lecture echouee pour <co2ice>"
-!         PRINT*, NF_STRERROR(ierr)
-!         CALL abort
-!      ENDIF
-c
+
       ierr = NF_INQ_VARID (nid, "emis", nvarid)
       IF (ierr .NE. NF_NOERR) THEN
@@ -718,81 +678,5 @@
          CALL abort
       ENDIF
-c
-cc Slab ocean
-      if(ok_slab_ocean) then
-      start=(/1,1,1,memo/)
-      count=(/imold+1,jmold+1,noceanmx,1/)
-
-       ierr=NF_INQ_VARID(nid,"tslab",nvarid)
-       IF (ierr.ne.NF_NOERR) then
-          PRINT*,"lect_start_archive: Cannot find <tslab>"
-       ENDIF
-#ifdef NC_DOUBLE
-      ierr=NF_GET_VARA_DOUBLE(nid,nvarid,start,count,tslabold)
-#else
-      ierr=NF_GET_VARA_REAL(nid,nvarid,start,count,tslabold)
-#endif
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "lect_start_archive: Lecture echouee pour <tslab>"
-      ENDIF
-
-
-c
-      start=(/1,1,memo,0/)
-      count=(/imold+1,jmold+1,1,0/)
-
-      ierr = NF_INQ_VARID (nid, "rnat", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "lect_start_archive: Le champ <rnat> est absent"
-      ENDIF
-#ifdef NC_DOUBLE
-      ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,count,rnatold)
-#else
-      ierr = NF_GET_VARA_REAL(nid, nvarid,start,count,rnatold)
-#endif
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "lect_start_archive: Lecture echouee pour <rnat>"
-      ENDIF
-c
-      ierr = NF_INQ_VARID (nid, "pctsrf_sic", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "lect_start_archive: Le champ <pctsrf_sic> est absent"
-      ENDIF
-#ifdef NC_DOUBLE
-      ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,count,pctsrf_sicold)
-#else
-      ierr = NF_GET_VARA_REAL(nid, nvarid,start,count,pctsrf_sicold)
-#endif
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "lect_start_archive: Lecture echouee pour <pctsrf_sic>"
-      ENDIF
-c
-      ierr = NF_INQ_VARID (nid, "tsea_ice", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "lect_start_archive: Le champ <tsea_ice> est absent"
-      ENDIF
-#ifdef NC_DOUBLE
-      ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,count,tsea_iceold)
-#else
-      ierr = NF_GET_VARA_REAL(nid, nvarid,start,count,tsea_iceold)
-#endif
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "lect_start_archive: Lecture echouee pour <tsea_ice>"
-      ENDIF
-c
-      ierr = NF_INQ_VARID (nid, "sea_ice", nvarid)
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "lect_start_archive: Le champ <sea_ice> est absent"
-      ENDIF
-#ifdef NC_DOUBLE
-      ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,count,sea_iceold)
-#else
-      ierr = NF_GET_VARA_REAL(nid, nvarid,start,count,sea_iceold)
-#endif
-      IF (ierr .NE. NF_NOERR) THEN
-         PRINT*, "lect_start_archive: Lecture echouee pour <sea_ice>"
-      ENDIF
-  
-      ENDIF! ok_slab_ocean
+
 c
       write(*,*)"lect_start_archive: rlonuold:"
@@ -813,13 +697,4 @@
 
       DO iq=1,nqtot
-          txt=trim(tname(iq))//"_surf"
-          if (txt.eq."h2o_vap") then
-            ! There is no surface tracer for h2o_vap;
-            ! "h2o_ice" should be loaded instead
-            txt="h2o_ice_surf"
-            write(*,*) 'lect_start_archive: loading surface tracer',
-     &                     ' h2o_ice instead of h2o_vap'
-          endif
-
         
         write(*,*) "lect_start_archive: loading tracer ",trim(txt)
@@ -1045,8 +920,4 @@
      &                   rlonuold,rlatvold,rlonu,rlatv)
 
-! CO2 ice is now in qsurf(igcm_co2_ice)
-!      call interp_horiz (co2iceold,co2ices,imold,jmold,iim,jjm,1,
-!     &                   rlonuold,rlatvold,rlonu,rlatv)
-
 c Temperature de surface
       call interp_horiz (tsurfold,tsurfs,imold,jmold,iim,jjm,1,
@@ -1088,20 +959,5 @@
          END DO
       END DO
-      co2icetotal = 0.
-      if (igcm_co2_ice.ne.0) then
-        ! recast surface CO2 ice on new grid
-        call interp_horiz(qsurfold(1,1,igcm_co2_ice),
-     &                  qsurfs(1,1,igcm_co2_ice),
-     &                  imold,jmold,iim,jjm,1,
-     &                  rlonuold,rlatvold,rlonu,rlatv)
-       DO j=1,jjp1
-         DO i=1,iim
-           !co2icetotal = co2icetotal + co2iceS(i,j)*aire(i,j)
-           co2icetotal=co2icetotal+qsurfS(i,j,igcm_co2_ice)*aire(i,j)
-         END DO
-       END DO
-      else
-        write(*,*) "Warning: No co2_ice surface tracer"
-      endif
+      
 
       write(*,*)
@@ -1110,11 +966,4 @@
       write (*,*) 'Ratio new atm./ old atm =', ptotal/ptotalold 
       write(*,*)
-      write(*,*)'Ancienne grille: masse de la glace CO2:',co2icetotalold
-      write(*,*)'Nouvelle grille: masse de la glace CO2:',co2icetotal
-      if (co2icetotalold.ne.0.) then
-      write(*,*)'Ratio new ice./old ice =',co2icetotal/co2icetotalold
-      endif
-      write(*,*)
-
 
       DO j=1,jjp1
@@ -1124,13 +973,4 @@
       END DO
 
-      if ( co2icetotalold.gt.0.) then 
-!         DO j=1,jjp1
-!            DO i=1,iip1
-!               co2iceS(i,j)=co2iceS(i,j) * co2icetotalold/co2icetotal
-!            END DO
-!         END DO
-        write(*,*) "Not enforcing conservation of surface CO2 ice"
-        write(*,*) " (should be OK when CO2 ice is a tracer)"
-      end if
 
 c-----------------------------------------------------------------------
@@ -1290,28 +1130,5 @@
 
 c-----------------------------------------------------------------------
-c 6.3   Slab Ocean :
-c-----------------------------------------------------------------------
-      call interp_horiz (tslabold,tslabs,imold,jmold,iim,jjm,noceanmx,
-     &                   rlonuold,rlatvold,rlonu,rlatv)
-      call gr_dyn_fi (noceanmx,iim+1,jjm+1,ngrid,tslabs,tslab)
-
-      call interp_horiz (rnatold,rnats,imold,jmold,iim,jjm,1,
-     &                   rlonuold,rlatvold,rlonu,rlatv)
-      call gr_dyn_fi (1,iim+1,jjm+1,ngrid,rnats,rnat)
-
-      call interp_horiz (pctsrf_sicold,pctsrf_sics,imold,jmold,iim,
-     &                   jjm,1,rlonuold,rlatvold,rlonu,rlatv)
-      call gr_dyn_fi (1,iim+1,jjm+1,ngrid,pctsrf_sics,pctsrf_sic)
-
-      call interp_horiz (tsea_iceold,tsea_ices,imold,jmold,iim,jjm,1,
-     &                   rlonuold,rlatvold,rlonu,rlatv)
-      call gr_dyn_fi (1,iim+1,jjm+1,ngrid,tsea_ices,tsea_ice)
-
-      call interp_horiz (sea_iceold,sea_ices,imold,jmold,iim,jjm,1,
-     &                   rlonuold,rlatvold,rlonu,rlatv)
-      call gr_dyn_fi (1,iim+1,jjm+1,ngrid,sea_ices,sea_ice)
-
-c-----------------------------------------------------------------------
-c 6.4 Variable 3d :
+c 6.3 Variable 3d :
 c-----------------------------------------------------------------------
       
@@ -1434,7 +1251,4 @@
          end do
       enddo
-      
-!      call gr_dyn_fi (1,iim+1,jjm+1,ngrid,co2ices,co2ice)
-! no need to transfer "co2ice" any more; it is in qsurf(igcm_co2_ice)
 
       endif !! if nqtot .ne. 0
@@ -1470,5 +1284,4 @@
       deallocate(phisold)
       deallocate(qold)
-      deallocate(co2iceold)
       deallocate(tsurfold)
       deallocate(emisold)
@@ -1482,9 +1295,4 @@
       deallocate(qsurfold)
       deallocate(var,varp1)
-      deallocate(tslabold)
-      deallocate(rnatold)
-      deallocate(pctsrf_sicold)
-      deallocate(tsea_iceold)
-      deallocate(sea_iceold)
 
 !      write(*,*)'lect_start_archive: END'
Index: trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/newstart.F
===================================================================
--- trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/newstart.F	(revision 1644)
+++ trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/newstart.F	(revision 1647)
@@ -19,5 +19,4 @@
      &                              is_master
       use infotrac, only: infotrac_init, nqtot, tname
-      USE tracer_h, ONLY: igcm_co2_ice, igcm_h2o_vap, igcm_h2o_ice
       USE comsoil_h, ONLY: nsoilmx, layer, mlayer, inertiedat
       USE surfdat_h, ONLY: phisfi, albedodat,
@@ -28,5 +27,4 @@
       use phyredem, only: physdem0, physdem1
       use iostart, only: open_startphy
-      use slab_ice_h, only:noceanmx
       use filtreg_mod, only: inifilr
       USE mod_const_mpi, ONLY: COMM_LMDZ
@@ -95,5 +93,4 @@
       REAL tsurf(ngridmx)        ! surface temperature
       REAL,ALLOCATABLE :: tsoil(:,:) ! soil temperature
-!      REAL co2ice(ngridmx)        ! CO2 ice layer
       REAL emis(ngridmx)        ! surface emissivity
       real emisread             ! added by RW
@@ -109,11 +106,5 @@
       real mugaz ! molar mass of the atmosphere
 
-      integer ierr 
-
-      REAL rnat(ngridmx)
-      REAL,ALLOCATABLE :: tslab(:,:) ! slab ocean temperature
-      REAL pctsrf_sic(ngridmx) ! sea ice cover
-      REAL tsea_ice(ngridmx) ! temperature sea_ice
-      REAL sea_ice(ngridmx) ! mass of sea ice (kg/m2)
+      integer ierr
 
 c Variables on the new grid along scalar points 
@@ -153,5 +144,4 @@
       logical :: flagtset=.false. ,  flagps0=.false.
       real val, val2, val3, val4 ! to store temporary variables
-      real :: iceith=2000 ! thermal inertia of subterranean ice
 
       INTEGER :: itau
@@ -165,9 +155,4 @@
 !     added by BC for equilibrium temperature startup
       real teque
-
-!     added by BC for cloud fraction setup
-      REAL hice(ngridmx),cloudfrac(ngridmx,llm)
-      REAL totalfrac(ngridmx)
-
 
 !     added by RW for nuketharsis
@@ -209,5 +194,4 @@
       allocate(tsoil(ngridmx,nsoilmx))
       allocate(ith(iip1,jjp1,nsoilmx),ithfi(ngridmx,nsoilmx))
-      allocate(tslab(ngridmx,nsoilmx))
       
 c=======================================================================
@@ -355,7 +339,5 @@
         CALL phyetat0 (ngridmx,llm,fichnom,tab0,Lmodif,nsoilmx,
      .        nqtot,day_ini,time,
-     .        tsurf,tsoil,emis,q2,qsurf,   !) ! temporary modif by RDW
-     .        cloudfrac,totalfrac,hice,rnat,pctsrf_sic,tslab,tsea_ice,
-     .        sea_ice)
+     .        tsurf,tsoil,emis,q2,qsurf)
 
         ! copy albedo and soil thermal inertia on (local) physics grid
@@ -549,6 +531,5 @@
      &   date,tsurf,tsoil,emis,q2,
      &   t,ucov,vcov,ps,teta,phisold_newgrid,q,qsurf,
-     &   surfith,nid,
-     &   rnat,pctsrf_sic,tslab,tsea_ice,sea_ice)
+     &   surfith,nid)
         write(*,*) "OK, read start_archive file"
         ! copy soil thermal inertia
@@ -784,9 +765,4 @@
         else if (trim(modif).eq.'ptot') then
 
-          ! check if we have a co2_ice surface tracer:
-          if (igcm_co2_ice.eq.0) then
-            write(*,*) " No surface CO2 ice !"
-            write(*,*) " only atmospheric pressure will be considered!"
-          endif
 c         calcul de la pression totale glace + atm actuelle
           patm=0.
@@ -800,13 +776,9 @@
                 patm = patm + ps(i,j)*aire(i,j)
                 airetot= airetot + aire(i,j)
-                if (igcm_co2_ice.ne.0) then
-                  !pcap = pcap + aire(i,j)*co2ice(ig)*g
-                  pcap = pcap + aire(i,j)*qsurf(ig,igcm_co2_ice)*g
-                endif
              ENDDO
           ENDDO
           ptoto = pcap + patm
 
-          print*,'Current total pressure at surface (co2 ice + atm) ',
+          print*,'Current total pressure at surface (atm) ',
      &       ptoto/airetot
 
@@ -1048,36 +1020,15 @@
        else if (trim(modif) .eq. 'wetstart') then
         ! check that there is indeed a water vapor tracer
-        if (igcm_h2o_vap.eq.0) then
+
           write(*,*) "No water vapour tracer! Can't use this option"
           stop
-        endif
-          DO l=1,llm
-            DO j=1,jjp1
-              DO i=1,iip1
-                q(i,j,l,igcm_h2o_vap)=150.e-6 * (rlatu(j)+pi/2.) / pi
-              ENDDO
-            ENDDO
-          ENDDO
-
-         write(*,*) 'Water mass mixing ratio at north pole='
-     *               ,q(1,1,1,igcm_h2o_vap)
-         write(*,*) '---------------------------south pole='
-     *               ,q(1,jjp1,1,igcm_h2o_vap)
 
 c      noglacier : remove tropical water ice (to initialize high res sim)
 c      --------------------------------------------------
         else if (trim(modif) .eq. 'noglacier') then
-           if (igcm_h2o_ice.eq.0) then
+        
              write(*,*) "No water ice tracer! Can't use this option"
              stop
-           endif
-           do ig=1,ngridmx
-             j=(ig-2)/iim +2
-              if(ig.eq.1) j=1
-              write(*,*) 'OK: remove surface ice for |lat|<45'
-              if (abs(rlatu(j)*180./pi).lt.45.) then
-                   qsurf(ig,igcm_h2o_ice)=0.
-              end if
-           end do
+
 
 
@@ -1085,267 +1036,56 @@
 c      --------------------------------------------------
         else if (trim(modif) .eq. 'watercapn') then
-           if (igcm_h2o_ice.eq.0) then
+        
              write(*,*) "No water ice tracer! Can't use this option"
              stop
-           endif
-
-           DO j=1,jjp1        
-              DO i=1,iim
-                 ig=1+(j-2)*iim +i
-                 if(j.eq.1) ig=1
-                 if(j.eq.jjp1) ig=ngridmx
-
-                 if (rlatu(j)*180./pi.gt.80.) then
-                    qsurf(ig,igcm_h2o_ice)=3.4e3
-                    !do isoil=1,nsoilmx
-                    !   ith(i,j,isoil) = 18000. ! thermal inertia
-                    !enddo
-                   write(*,*)'     ==> Ice mesh North boundary (deg)= ',
-     &                   rlatv(j-1)*180./pi
-                 end if
-              ENDDO
-           ENDDO
-           CALL gr_dyn_fi(nsoilmx,iip1,jjp1,ngridmx,ith,ithfi)
-
-c$$$           do ig=1,ngridmx
-c$$$             j=(ig-2)/iim +2
-c$$$              if(ig.eq.1) j=1
-c$$$              if (rlatu(j)*180./pi.gt.80.) then
-c$$$
-c$$$                   qsurf(ig,igcm_h2o_ice)=1.e5
-c$$$                   qsurf(ig,igcm_h2o_vap)=0.0!1.e5
-c$$$
-c$$$                   write(*,*) 'ig=',ig,'    H2O ice mass (kg/m2)= ',
-c$$$     &              qsurf(ig,igcm_h2o_ice)
-c$$$
-c$$$                   write(*,*)'     ==> Ice mesh South boundary (deg)= ',
-c$$$     &              rlatv(j)*180./pi
-c$$$                end if
-c$$$           enddo
 
 c      watercaps : H20 ice on permanent southern cap
 c      -------------------------------------------------
         else if (trim(modif) .eq. 'watercaps') then
-           if (igcm_h2o_ice.eq.0) then
+
               write(*,*) "No water ice tracer! Can't use this option"
               stop
-           endif
-
-           DO j=1,jjp1        
-              DO i=1,iim
-                 ig=1+(j-2)*iim +i
-                 if(j.eq.1) ig=1
-                 if(j.eq.jjp1) ig=ngridmx
-
-                 if (rlatu(j)*180./pi.lt.-80.) then
-                    qsurf(ig,igcm_h2o_ice)=3.4e3
-                    !do isoil=1,nsoilmx
-                    !   ith(i,j,isoil) = 18000. ! thermal inertia
-                    !enddo
-                   write(*,*)'     ==> Ice mesh South boundary (deg)= ',
-     &                   rlatv(j-1)*180./pi
-                 end if
-              ENDDO
-           ENDDO
-           CALL gr_dyn_fi(nsoilmx,iip1,jjp1,ngridmx,ith,ithfi)
-
-c$$$           do ig=1,ngridmx
-c$$$               j=(ig-2)/iim +2
-c$$$               if(ig.eq.1) j=1
-c$$$               if (rlatu(j)*180./pi.lt.-80.) then
-c$$$                  qsurf(ig,igcm_h2o_ice)=1.e5
-c$$$                  qsurf(ig,igcm_h2o_vap)=0.0 !1.e5
-c$$$
-c$$$                  write(*,*) 'ig=',ig,'   H2O ice mass (kg/m2)= ',
-c$$$     &                 qsurf(ig,igcm_h2o_ice)
-c$$$                  write(*,*)'     ==> Ice mesh North boundary (deg)= ',
-c$$$     &                 rlatv(j-1)*180./pi
-c$$$               end if
-c$$$           enddo
-
 
 c       noacglac : H2O ice across highest terrain
 c       --------------------------------------------
         else if (trim(modif) .eq. 'noacglac') then
-           if (igcm_h2o_ice.eq.0) then
+
              write(*,*) "No water ice tracer! Can't use this option"
              stop
-           endif
-          DO j=1,jjp1        
-             DO i=1,iim
-                ig=1+(j-2)*iim +i
-                if(j.eq.1) ig=1
-                if(j.eq.jjp1) ig=ngridmx
-
-                if(phis(i,j).gt.1000.*g)then
-                    alb(i,j) = 0.5 ! snow value
-                    do isoil=1,nsoilmx
-                       ith(i,j,isoil) = 18000. ! thermal inertia
-                       ! this leads to rnat set to 'ocean' in physiq.F90
-                       ! actually no, because it is soil not surface
-                    enddo
-                endif
-             ENDDO
-          ENDDO
-          CALL gr_dyn_fi(nsoilmx,iip1,jjp1,ngridmx,ith,ithfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,alb,albfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,phis,phisfi)
-
-
 
 c       oborealis : H2O oceans across Vastitas Borealis
 c       -----------------------------------------------
         else if (trim(modif) .eq. 'oborealis') then
-           if (igcm_h2o_ice.eq.0) then
+
              write(*,*) "No water ice tracer! Can't use this option"
              stop
-           endif
-          DO j=1,jjp1        
-             DO i=1,iim
-                ig=1+(j-2)*iim +i
-                if(j.eq.1) ig=1
-                if(j.eq.jjp1) ig=ngridmx
-
-                if(phis(i,j).lt.-4000.*g)then
-!                if( (phis(i,j).lt.-4000.*g)
-!     &               .and. (rlatu(j)*180./pi.lt.0.) )then ! south hemisphere only
-!                    phis(i,j)=-8000.0*g ! proper ocean
-                    phis(i,j)=-4000.0*g ! scanty ocean
-
-                    alb(i,j) = 0.07 ! oceanic value
-                    do isoil=1,nsoilmx
-                       ith(i,j,isoil) = 18000. ! thermal inertia
-                       ! this leads to rnat set to 'ocean' in physiq.F90
-                       ! actually no, because it is soil not surface
-                    enddo
-                endif
-             ENDDO
-          ENDDO
-          CALL gr_dyn_fi(nsoilmx,iip1,jjp1,ngridmx,ith,ithfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,alb,albfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,phis,phisfi)
-
+             
 c       iborealis : H2O ice in Northern plains
 c       --------------------------------------
         else if (trim(modif) .eq. 'iborealis') then
-           if (igcm_h2o_ice.eq.0) then
+
              write(*,*) "No water ice tracer! Can't use this option"
              stop
-           endif
-          DO j=1,jjp1        
-             DO i=1,iim
-                ig=1+(j-2)*iim +i
-                if(j.eq.1) ig=1
-                if(j.eq.jjp1) ig=ngridmx
-
-                if(phis(i,j).lt.-4000.*g)then
-                   !qsurf(ig,igcm_h2o_ice)=1.e3
-                   qsurf(ig,igcm_h2o_ice)=241.4 ! to make total 33 kg m^-2
-                endif
-             ENDDO
-          ENDDO
-          CALL gr_dyn_fi(nsoilmx,iip1,jjp1,ngridmx,ith,ithfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,alb,albfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,phis,phisfi)
-
 
 c       oceanball : H2O liquid everywhere
 c       ----------------------------
         else if (trim(modif) .eq. 'oceanball') then
-           if (igcm_h2o_ice.eq.0) then
+
              write(*,*) "No water ice tracer! Can't use this option"
              stop
-           endif
-          DO j=1,jjp1        
-             DO i=1,iim
-                ig=1+(j-2)*iim +i
-                if(j.eq.1) ig=1
-                if(j.eq.jjp1) ig=ngridmx
-
-                qsurf(ig,igcm_h2o_vap)=0.0    ! for ocean, this is infinite source
-                qsurf(ig,igcm_h2o_ice)=0.0
-                alb(i,j) = 0.07 ! ocean value
-
-                do isoil=1,nsoilmx
-                   ith(i,j,isoil) = 18000. ! thermal inertia
-                   !ith(i,j,isoil) = 50. ! extremely low for test
-                   ! this leads to rnat set to 'ocean' in physiq.F90
-                enddo
-
-             ENDDO
-          ENDDO
-          CALL gr_dyn_fi(nsoilmx,iip1,jjp1,ngridmx,ith,ithfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,alb,albfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,phis,phisfi)
 
 c       iceball : H2O ice everywhere
 c       ----------------------------
         else if (trim(modif) .eq. 'iceball') then
-           if (igcm_h2o_ice.eq.0) then
+
              write(*,*) "No water ice tracer! Can't use this option"
              stop
-           endif
-          DO j=1,jjp1        
-             DO i=1,iim
-                ig=1+(j-2)*iim +i
-                if(j.eq.1) ig=1
-                if(j.eq.jjp1) ig=ngridmx
-
-                qsurf(ig,igcm_h2o_vap)=-50.    ! for ocean, this is infinite source
-                qsurf(ig,igcm_h2o_ice)=50.     ! == to 0.05 m of oceanic ice
-                alb(i,j) = 0.6 ! ice albedo value
-
-                do isoil=1,nsoilmx
-                   !ith(i,j,isoil) = 18000. ! thermal inertia
-                   ! this leads to rnat set to 'ocean' in physiq.F90
-                enddo
-
-             ENDDO
-          ENDDO
-          CALL gr_dyn_fi(nsoilmx,iip1,jjp1,ngridmx,ith,ithfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,alb,albfi)
 
 c       supercontinent : H2O ice everywhere
 c       ----------------------------
         else if (trim(modif) .eq. 'supercontinent') then
-          write(*,*) 'Minimum longitude (-180,180)'
-          read(*,*) val
-          write(*,*) 'Maximum longitude (-180,180)'
-          read(*,*) val2
-          write(*,*) 'Minimum latitude (-90,90)'
-          read(*,*) val3
-          write(*,*) 'Maximum latitude (-90,90)'
-          read(*,*) val4
-
-          do j=1,jjp1
-            do i=1,iip1
-              ig=1+(j-2)*iim +i
-              if(j.eq.1) ig=1
-              if(j.eq.jjp1) ig=ngridmx
-
-c             Supercontinent:
-              if (((rlatu(j)*180./pi.le.val4).and.
-     &            (rlatu(j)*180./pi.ge.val3).and.
-     &            (rlonv(i)*180./pi.le.val2).and.
-     &            (rlonv(i)*180./pi.ge.val))) then
-
-                rnat(ig)=1.
-                alb(i,j) = 0.3
-                do isoil=1,nsoilmx
-                  ith(i,j,isoil) = 2000.
-                enddo
-c             Ocean:
-              else
-                rnat(ig)=0.
-                alb(i,j) = 0.07
-                do isoil=1,nsoilmx
-                  ith(i,j,isoil) = 18000.
-                enddo
-              end if
-
-            enddo
-          enddo
-          CALL gr_dyn_fi(nsoilmx,iip1,jjp1,ngridmx,ith,ithfi)
-          CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,alb,albfi)
+ 
+             write(*,*) "No water ice tracer! Can't use this option"
+             stop
 
 c       isotherm : Isothermal temperatures and no winds
@@ -1445,14 +1185,6 @@
 c       ------------------------------------------------
         else if (trim(modif) .eq. 'co2ice=0') then
-         ! check that there is indeed a co2_ice tracer ...
-          if (igcm_co2_ice.ne.0) then
-           do ig=1,ngridmx
-              !co2ice(ig)=0
-              qsurf(ig,igcm_co2_ice)=0
-              emis(ig)=emis(ngridmx/2)
-           end do
-          else
             write(*,*) "Can't remove CO2 ice!! (no co2_ice tracer)"
-          endif
+        
         
 !       therm_ini_s: (re)-set soil thermal inertia to reference surface values
@@ -1484,26 +1216,7 @@
 c       ------------------------------------------------
         else if (trim(modif) .eq. 'slab_ocean_0') then
-        write(*,*)'OK: initialisation of slab ocean' 
-
-      DO ig=1, ngridmx
-         rnat(ig)=1.
-         tslab(ig,1)=0.
-         tslab(ig,2)=0.
-         tsea_ice(ig)=0.
-         sea_ice(ig)=0.
-         pctsrf_sic(ig)=0.
-         
-         if(ithfi(ig,1).GT.10000.)then
-           rnat(ig)=0.
-           phisfi(ig)=0.
-           tsea_ice(ig)=273.15-1.8
-           tslab(ig,1)=tsurf(ig)
-           tslab(ig,2)=tsurf(ig)!*2./3.+(273.15-1.8)/3.
-          endif
-
-      ENDDO
-          CALL gr_fi_dyn(1,ngridmx,iip1,jjp1,phisfi,phis)
-
-
+        
+           write(*,*) "No ocean yet on Titan ! Can't use this option"
+           stop
 
         else
@@ -1517,28 +1230,5 @@
  999  continue
 
- 
-c=======================================================================
-c   Initialisation for cloud fraction and oceanic ice (added by BC 2010)
-c=======================================================================
-      DO ig=1, ngridmx
-         totalfrac(ig)=0.5
-         DO l=1,llm
-            cloudfrac(ig,l)=0.5
-         ENDDO
-! Ehouarn, march 2012: also add some initialisation for hice
-         hice(ig)=0.0
-      ENDDO
       
-c========================================================
-
-!      DO ig=1,ngridmx
-!         IF(tsurf(ig) .LT. 273.16-1.8) THEN
-!            hice(ig)=(273.16-1.8-tsurf(ig))/(273.16-1.8-240)*1.
-!            hice(ig)=min(hice(ig),1.0)
-!         ENDIF
-!      ENDDO
-
-
-
 
 c=======================================================================
@@ -1647,7 +1337,5 @@
       call physdem1("restartfi.nc",nsoilmx,ngridmx,llm,nqtot,
      &                dtphys,real(day_ini),
-     &                tsurf,tsoil,emis,q2,qsurf,
-     &                cloudfrac,totalfrac,hice,
-     &                rnat,pctsrf_sic,tslab,tsea_ice,sea_ice)
+     &                tsurf,tsoil,emis,q2,qsurf)
 
 c=======================================================================
Index: trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/start2archive.F
===================================================================
--- trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/start2archive.F	(revision 1644)
+++ trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/start2archive.F	(revision 1647)
@@ -25,5 +25,4 @@
 !      use control_mod
 !      use comgeomphy, only: initcomgeomphy
-      use slab_ice_h, only: noceanmx
 ! to use  'getin'
       USE ioipsl_getincom
@@ -31,5 +30,4 @@
       USE mod_const_mpi, ONLY: COMM_LMDZ
       USE control_mod, only: planet_type
-      USE callkeys_mod, ONLY: ok_slab_ocean
       use filtreg_mod, only: inifilr
       USE comvert_mod, ONLY: ap,bp
@@ -72,5 +70,4 @@
       REAL tsurf(ngridmx)        ! Surface temperature
       REAL,ALLOCATABLE :: tsoil(:,:) ! Soil temperature
-      REAL co2ice(ngridmx)        ! CO2 ice layer
       REAL q2(ngridmx,llm+1)
       REAL,ALLOCATABLE :: qsurf(:,:)
@@ -82,12 +79,4 @@
       INTEGER*4 day_ini_fi
 
-!     added by FF for cloud fraction setup
-      REAL hice(ngridmx)
-      REAL cloudfrac(ngridmx,llm),totalcloudfrac(ngridmx)
-
-!     added by BC for slab ocean
-      REAL rnat(ngridmx),pctsrf_sic(ngridmx),sea_ice(ngridmx)
-      REAL tslab(ngridmx,noceanmx),tsea_ice(ngridmx)
-
 
 c Variable naturelle / grille scalaire
@@ -97,16 +86,7 @@
       REAL,ALLOCATABLE :: tsoilS(:,:)
       REAL,ALLOCATABLE :: ithS(:,:) ! Soil Thermal Inertia
-      REAL co2iceS(ip1jmp1)
       REAL q2S(ip1jmp1,llm+1)
       REAL,ALLOCATABLE :: qsurfS(:,:)
       REAL emisS(ip1jmp1)
-
-!     added by FF for cloud fraction setup
-      REAL hiceS(ip1jmp1)
-      REAL cloudfracS(ip1jmp1,llm),totalcloudfracS(ip1jmp1)
-
-!     added by BC for slab ocean
-      REAL rnatS(ip1jmp1),pctsrf_sicS(ip1jmp1),sea_iceS(ip1jmp1)
-      REAL tslabS(ip1jmp1,noceanmx),tsea_iceS(ip1jmp1)
 
 
@@ -120,5 +100,5 @@
       INTEGER Lmodif
 
-      REAL ptotal, co2icetotal
+      REAL ptotal
       REAL timedyn,timefi !fraction du jour dans start, startfi
       REAL date
@@ -236,10 +216,5 @@
       CALL phyetat0 (ngridmx,llm,fichnom,0,Lmodif,nsoilmx,nqtot,
      .      day_ini_fi,timefi,
-     .      tsurf,tsoil,emis,q2,qsurf,
-!       change FF 05/2011
-     .       cloudfrac,totalcloudfrac,hice,
-!       change BC 05/2014
-     .       rnat,pctsrf_sic,tslab,tsea_ice,sea_ice)
-
+     .      tsurf,tsoil,emis,q2,qsurf)
 
 
@@ -328,5 +303,4 @@
 c
 c tsurf --> tsurfS
-c co2ice --> co2iceS
 c tsoil --> tsoilS
 c emis --> emisS
@@ -337,5 +311,4 @@
 
       call gr_fi_dyn(1,ngridmx,iip1,jjp1,tsurf,tsurfS)
-!      call gr_fi_dyn(1,ngridmx,iip1,jjp1,co2ice,co2iceS)
       call gr_fi_dyn(nsoilmx,ngridmx,iip1,jjp1,tsoil,tsoilS)
       ! Note: thermal inertia "inertiedat" is in comsoil.h
@@ -344,13 +317,4 @@
       call gr_fi_dyn(llm+1,ngridmx,iip1,jjp1,q2,q2S)
       call gr_fi_dyn(nqtot,ngridmx,iip1,jjp1,qsurf,qsurfS)
-      call gr_fi_dyn(llm,ngridmx,iip1,jjp1,cloudfrac,cloudfracS)
-      call gr_fi_dyn(1,ngridmx,iip1,jjp1,hice,hiceS)
-      call gr_fi_dyn(1,ngridmx,iip1,jjp1,totalcloudfrac,totalcloudfracS)
-
-      call gr_fi_dyn(1,ngridmx,iip1,jjp1,rnat,rnatS)
-      call gr_fi_dyn(1,ngridmx,iip1,jjp1,pctsrf_sic,pctsrf_sicS)
-      call gr_fi_dyn(1,ngridmx,iip1,jjp1,tsea_ice,tsea_iceS)
-      call gr_fi_dyn(1,ngridmx,iip1,jjp1,sea_ice,sea_iceS)
-      call gr_fi_dyn(noceanmx,ngridmx,iip1,jjp1,tslab,tslabS)
 
 c=======================================================================
@@ -359,21 +323,17 @@
 
       ptotal =  0.
-      co2icetotal = 0.
       DO j=1,jjp1
          DO i=1,iim
            ptotal=ptotal+aire(i+(iim+1)*(j-1))*ps(i+(iim+1)*(j-1))/g
-!           co2icetotal = co2icetotal + 
-!     &            co2iceS(i+(iim+1)*(j-1))*aire(i+(iim+1)*(j-1))
          ENDDO
       ENDDO
       write(*,*)'Ancienne grille : masse de l''atm :',ptotal
-!      write(*,*)'Ancienne grille : masse de la glace CO2 :',co2icetotal
-
-c-----------------------------------------------------------------------
-c Passage de "ptotal" et "co2icetotal" par tab_cntrl_fi
+
+c-----------------------------------------------------------------------
+c Passage de "ptotal"  par tab_cntrl_fi
 c-----------------------------------------------------------------------
 
       tab_cntrl_fi(49) = ptotal
-      tab_cntrl_fi(50) = co2icetotal
+      tab_cntrl_fi(50) = 0.
 
 c=======================================================================
@@ -430,5 +390,5 @@
 
 c-----------------------------------------------------------------------
-c Ecriture des champs  (co2ice,emis,ps,Tsurf,T,u,v,q2,q,qsurf)
+c Ecriture des champs  (emis,ps,Tsurf,T,u,v,q2,q,qsurf)
 c-----------------------------------------------------------------------
 c ATTENTION: q2 a une couche de plus!!!!
@@ -439,6 +399,4 @@
 c-----------------------------------------------------------------------
 
-!      call write_archive(nid,ntime,'co2ice','couche de glace co2',
-!     &  'kg/m2',2,co2iceS)
       call write_archive(nid,ntime,'emis','grd emis',' ',2,emisS)
       call write_archive(nid,ntime,'ps','Psurf','Pa',2,ps)
@@ -497,43 +455,5 @@
 ! Note: no need to write volcapa, it is stored in "controle" table
 
-c-----------------------------------------------------------------------
-c Ecriture du champs  cloudfrac,hice,totalcloudfrac
-c-----------------------------------------------------------------------
-      call write_archive(nid,ntime,'hice',
-     &         'Height of oceanic ice','m',2,hiceS)
-      call write_archive(nid,ntime,'totalcloudfrac',
-     &        'Total cloud Fraction','',2,totalcloudfracS)
-      call write_archive(nid,ntime,'cloudfrac'
-     &        ,'Cloud fraction','',3,cloudfracS)
-
-c-----------------------------------------------------------------------
-c Slab ocean
-c-----------------------------------------------------------------------
-      OPEN(99,file='callphys.def',status='old',form='formatted'
-     &     ,iostat=ierr)
-      CLOSE(99)
-
-      IF(ierr.EQ.0) THEN
-
-
-         write(*,*) "Use slab-ocean ?"
-         ok_slab_ocean=.false.         ! default value
-         call getin("ok_slab_ocean",ok_slab_ocean)
-         write(*,*) "ok_slab_ocean = ",ok_slab_ocean
-
-      if(ok_slab_ocean) then
-      call write_archive(nid,ntime,'rnat'
-     &        ,'rnat','',2,rnatS)
-      call write_archive(nid,ntime,'pctsrf_sic'
-     &        ,'pctsrf_sic','',2,pctsrf_sicS)
-      call write_archive(nid,ntime,'sea_ice'
-     &        ,'sea_ice','',2,sea_iceS)
-      call write_archive(nid,ntime,'tslab'
-     &        ,'tslab','',-2,tslabS)
-      call write_archive(nid,ntime,'tsea_ice'
-     &        ,'tsea_ice','',2,tsea_iceS)
-      endif !ok_slab_ocean
-      ENDIF
-c-----------------------------------------------------------------------
+
 c Fin 
 c-----------------------------------------------------------------------
Index: trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/write_archive.F
===================================================================
--- trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/write_archive.F	(revision 1644)
+++ trunk/LMDZ.TITAN/libf/dynphy_lonlat/phytitan/write_archive.F	(revision 1647)
@@ -33,5 +33,4 @@
 
       use comsoil_h, only: nsoilmx
-      use slab_ice_h, only: noceanmx
 
       implicit none
@@ -195,5 +194,5 @@
         edges(1)=iip1
         edges(2)=jjp1
-        edges(3)=noceanmx
+        edges(3)=1 ! JVO2017 : was noceanmx before -> set to 1
         edges(4)=1
 #ifdef NC_DOUBLE
