Index: /trunk/libf/phytitan/YOMCST.h
===================================================================
--- /trunk/libf/phytitan/YOMCST.h	(revision 101)
+++ /trunk/libf/phytitan/YOMCST.h	(revision 102)
@@ -10,9 +10,9 @@
       REAL RDAY,REA,REPSM,RSIYEA,RSIDAY,ROMEGA
 ! A1.1.bis Constantes concernant l'orbite de la Terre:
-      REAL :: R_ecc, R_peri, R_incl
+      REAL R_ecc, R_peri, R_incl
 ! A1.2 Geoide
       REAL RA,RG,R1SA
 ! A1.3 Radiation
-!     REAL :: RSIGMA,RI0
+!     REAL RSIGMA,RI0
       REAL RSIGMA
 ! A1.4 Thermodynamic gas phase
@@ -31,5 +31,5 @@
 !
       COMMON/YOMCST/RPI ,RCLUM, RHPLA, RKBOL, RNAVO ,RDAY  ,REA         &
-     & ,REPSM ,RSIYEA,RSIDAY,ROMEGA , R_ecc, R_peri, R_incl,            &
+     & ,REPSM ,RSIYEA,RSIDAY,ROMEGA , R_ecc, R_peri, R_incl             &
      & ,RA    ,RG ,R1SA                                                 &
      & ,RSIGMA,R ,RMD   ,RMV   ,RD    ,RV    ,RCPD ,RCPV,RCVD           &
Index: /trunk/libf/phytitan/aaam_bud.F
===================================================================
--- /trunk/libf/phytitan/aaam_bud.F	(revision 101)
+++ /trunk/libf/phytitan/aaam_bud.F	(revision 102)
@@ -2,8 +2,9 @@
      i                   rea,rg,ome,      
      i                   plat,plon,phis,
-     i                   dragu,liftu,phyu,
-     i                   dragv,liftv,phyv,
+     i                   dragu,liftu,clu,
+     i                   dragv,liftv,clv,
      i                   p, u, v)
 c
+      use dimphy
       implicit none
 c======================================================================
@@ -17,6 +18,6 @@
 C    C               ra,rg,romega,
 C    C               rlat,rlon,pphis,
-C    C               zustrdr,zustrli,zustrph,
-C    C               zvstrdr,zvstrli,zvstrph,
+C    C               zustrdr,zustrli,zustrcl,
+C    C               zvstrdr,zvstrli,zvstrcl,
 C    C               paprs,u,v)
 C
@@ -39,8 +40,8 @@
 c dragu---input-R-orodrag stress (zonal)
 c liftu---input-R-orolift stress (zonal)
-c phyu----input-R-Stress total de la physique (zonal)
+c clu-----input-R-Boundary layer stress (zonal)
 c dragv---input-R-orodrag stress (Meridional)
 c liftv---input-R-orolift stress (Meridional)
-c phyv----input-R-Stress total de la physique (Meridional)
+c clv-----input-R-Boundary layer stress (Meridional)
 c p-------input-R-Pressure (Pa) at model half levels
 c u-------input-R-Horizontal wind (m/s)
@@ -78,5 +79,4 @@
 
 #include "dimensions.h"
-#include "dimphy.h"
 c
 c ARGUMENTS
@@ -85,6 +85,6 @@
       REAL rjour,rsec,rea,rg,ome
       REAL plat(nlon),plon(nlon),phis(nlon)
-      REAL dragu(nlon),liftu(nlon),phyu(nlon)             
-      REAL dragv(nlon),liftv(nlon),phyv(nlon)             
+      REAL dragu(nlon),liftu(nlon),clu(nlon)             
+      REAL dragv(nlon),liftv(nlon),clv(nlon)             
       REAL p(nlon,nlev+1), u(nlon,nlev), v(nlon,nlev)
 c
@@ -148,6 +148,6 @@
           ssou(i,1)=dragu(l)+liftu(l)
           ssov(i,1)=dragv(l)+liftv(l)
-          blsu(i,1)=phyu(l)-dragu(l)-liftu(l)
-          blsv(i,1)=phyv(l)-dragv(l)-liftv(l)
+          blsu(i,1)=clu(l)
+          blsv(i,1)=clv(l)
 
         enddo
@@ -162,6 +162,6 @@
           ssou(iim+1,j)=dragu(l+1)+liftu(l+1)
           ssov(iim+1,j)=dragv(l+1)+liftv(l+1)
-          blsu(iim+1,j)=phyu(l+1)-dragu(l+1)-liftu(l+1)
-          blsv(iim+1,j)=phyv(l+1)-dragv(l+1)-liftv(l+1)
+          blsu(iim+1,j)=clu(l+1)
+          blsv(iim+1,j)=clv(l+1)
       zlon(iim+1)=-plon(l+1)*xpi/180.
       zlat(j)=plat(l+1)*xpi/180.
@@ -182,6 +182,6 @@
           ssou(i,j)=dragu(l)+liftu(l)
           ssov(i,j)=dragv(l)+liftv(l)
-          blsu(i,j)=phyu(l)-dragu(l)-liftu(l)
-          blsv(i,j)=phyv(l)-dragv(l)-liftv(l)
+          blsu(i,j)=clu(l)
+          blsv(i,j)=clv(l)
       zlon(i)=plon(l)*xpi/180.
 
@@ -214,6 +214,6 @@
           ssou(i,jjm+1)=dragu(l)+liftu(l)
           ssov(i,jjm+1)=dragv(l)+liftv(l)
-          blsu(i,jjm+1)=phyu(l)-dragu(l)-liftu(l)
-          blsv(i,jjm+1)=phyv(l)-dragv(l)-liftv(l)
+          blsu(i,jjm+1)=clu(l)
+          blsv(i,jjm+1)=clv(l)
       ub(i,jjm+1)=ub(1,jjm+1)                               
       vb(i,jjm+1)=vb(1,jjm+1)                                
@@ -333,5 +333,5 @@
      c      tmou(3)/hadley,tsso(3)/hadley,tbls(3)/hadley 
 100   format(F12.5,15(1x,F12.5))
-c100   format(F12.5,5(1x,F12.5))
+c00   format(F12.5,5(1x,F12.5))
 
       write(iam+1,*)((zs(i,j),i=1,iim),j=1,jjm+1)
Index: /trunk/libf/phytitan/ajsec.F
===================================================================
--- /trunk/libf/phytitan/ajsec.F	(revision 101)
+++ /trunk/libf/phytitan/ajsec.F	(revision 102)
@@ -3,6 +3,8 @@
 !
 ! ADAPTATION GCM POUR CP(T)
-      SUBROUTINE ajsec(paprs, pplay, ppk, tfi, ufi, vfi, nq, qfi,
+      SUBROUTINE ajsec(paprs, pplay, ppk, tfi, ufi, vfi, nq, qfi, 
      .                             d_tfi, d_ufi, d_vfi, d_qfi)
+
+      use dimphy
       IMPLICIT none
 c======================================================================
@@ -23,8 +25,6 @@
 c d_vfi-----output-R-Incrementation du vent meridien
 c d_qfi-----output-R-Incrementation des traceurs
-
 c======================================================================
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
       REAL paprs(klon,klev+1), pplay(klon,klev)
@@ -36,7 +36,5 @@
       REAL qfi(klon,klev,nq), d_qfi(klon,klev,nq)
 c
-      INTEGER limbas, limhau ! les couches a ajuster
-ccc      PARAMETER (limbas=klev-3, limhau=klev)
-      PARAMETER (limbas=1, limhau=klev)
+      INTEGER,save :: limbas, limhau ! les couches a ajuster
 c
       REAL zh(klon,klev)
@@ -51,4 +49,7 @@
 c Initialisation:
 c
+      limbas=1
+      limhau=klev
+
       DO k = 1, klev
       DO i = 1, klon
@@ -56,8 +57,8 @@
          d_ufi(i,k) = 0.0
          d_vfi(i,k) = 0.0
-	 d_qfi(i,k,:) = 0.0
+         d_qfi(i,k,:) = 0.0
          zu(i,k)    = ufi(i,k)
          zv(i,k)    = vfi(i,k)
-	 zq(i,k,:)  = qfi(i,k,:)
+         zq(i,k,:)  = qfi(i,k,:)
       ENDDO
       ENDDO
@@ -119,5 +120,5 @@
               zum=0.
               zvm=0.
-	      zqm=0.
+              zqm=0.
               DO k = k1, k2
                 zalpha=zalpha+ABS(zh(i,k)-hm)*zdp(i,k)
@@ -125,14 +126,14 @@
                 zum=zum+ufi(i,k)*zdp(i,k)
                 zvm=zvm+vfi(i,k)*zdp(i,k)
-		do iq=1,nq
-		   zqm(iq)=zqm(iq)+qfi(i,k,iq)*zdp(i,k)
-		enddo
+                do iq=1,nq
+                  zqm(iq)=zqm(iq)+qfi(i,k,iq)*zdp(i,k)
+                enddo
               ENDDO
               zalpha=zalpha/(hm*(paprs(i,k1)-paprs(i,k2+1)))
               zum=zum/(paprs(i,k1)-paprs(i,k2+1))
               zvm=zvm/(paprs(i,k1)-paprs(i,k2+1))
-	      do iq=1,nq
-	         zqm(iq)=zqm(iq)/(paprs(i,k1)-paprs(i,k2+1))
-	      enddo
+              do iq=1,nq
+                zqm(iq)=zqm(iq)/(paprs(i,k1)-paprs(i,k2+1))
+              enddo
 
               IF(zalpha.GT.1.) THEN
@@ -153,5 +154,5 @@
                  zv(i,k)=vfi(i,k)+zalpha*(zvm-vfi(i,k))
                  do iq=1,nq
-		   zq(i,k,iq)=qfi(i,k,iq)+zalpha*(zqm(iq)-qfi(i,k,iq))
+                   zq(i,k,iq)=qfi(i,k,iq)+zalpha*(zqm(iq)-qfi(i,k,iq))
                  enddo
               ENDDO
@@ -187,5 +188,5 @@
          do iq=1,nq
            d_qfi(i,k,iq) = 0.0
-         enddo 
+         enddo
       ENDDO
       ENDDO
@@ -200,5 +201,5 @@
          do iq=1,nq
            d_qfi(i,k,iq) = 0.0
-         enddo 
+         enddo
       ENDDO
       ENDDO
Index: /trunk/libf/phytitan/ballon.F
===================================================================
--- /trunk/libf/phytitan/ballon.F	(revision 101)
+++ /trunk/libf/phytitan/ballon.F	(revision 102)
@@ -2,4 +2,5 @@
      i                   temp, p, u, v, geop)
 
+      use dimphy
       implicit none
 
@@ -94,5 +95,4 @@
 
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 c
@@ -255,4 +255,5 @@
      .                       phi,lam,ubal,vbal,abal)
 
+      use dimphy
       implicit none
 
@@ -283,5 +284,4 @@
 
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 c
@@ -354,4 +354,5 @@
       subroutine pos_implem(phi,lam,ubal,vbal,dt)
 
+      use dimphy
       implicit none
 
@@ -370,5 +371,4 @@
 
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 c
Index: /trunk/libf/phytitan/calchim.F
===================================================================
--- /trunk/libf/phytitan/calchim.F	(revision 101)
+++ /trunk/libf/phytitan/calchim.F	(revision 102)
@@ -2,9 +2,8 @@
      .                   ctemp,cplay,cplev,
      .                   dqyc)
-      implicit none
       
 c-------------------------------------------------
 c  
-c     Introduction d'une routine chimique
+c     Introduction d une routine chimique
 c
 c     Auteur: S. Lebonnois,  01/2000 | 09/2003 
@@ -13,7 +12,7 @@
 c-------------------------------------------------
 c
-
+      use dimphy
+      implicit none
 #include "dimensions.h"
-#include "dimphy.h"
 #include "clesphys.h"
 #include "paramet.h"
Index: unk/libf/phytitan/calchim.dbleprec
===================================================================
--- /trunk/libf/phytitan/calchim.dbleprec	(revision 101)
+++ 	(revision )
@@ -1,565 +1,0 @@
-      SUBROUTINE calchim(ny,qy_c,nomqy_c,declin_rad,ls_rad,dtchim,
-     .                   ctemp,cplay,cplev,
-     .                   
-     .
-     .                   dqyc)
-      implicit none
-      
-c-------------------------------------------------
-c  
-c     Introduction d'une routine chimique
-c
-c     Auteur: S. Lebonnois,  01/2000 | 09/2003 
-c            adaptation pour Titan 3D: 02/2009
-c
-c-------------------------------------------------
-c
-
-#include "dimensions.h"
-#include "dimphy.h"
-#include "clesphys.h"
-#include "paramet.h"
-#include "YOMCST.h"
-
-#include "titan_for.h"  
-!!!  doit etre en accord avec titan.h
-#include "aerprod.h"
-
-c    Arguments
-c    ---------
-
-      INTEGER      ny                     ! nb de composes (nqmax-nmicro) 
-      REAL         qy_c(jjm+1,klev,NC)    ! Especes chimiques apres adv.+diss.
-      character*10 nomqy_c(NC+1)          ! Noms des especes chimiques
-      REAL         declin_rad,ls_rad      ! declinaison et long solaire en radians
-      REAL         dtchim                 ! pas de temps chimie
-      REAL         ctemp(jjm+1,klev)      ! Temperature
-      REAL         cplay(jjm+1,klev)      ! pression (Pa)
-      REAL         cplev(jjm+1,klev)      ! pression intercouches (Pa)
-      
-      REAL         dqyc(jjm+1,klev,nqmx)  ! Tendances especes chimiques (nqmx, mais en fait NC...)
-      
-c    Local variables :
-c    -----------------
-c variables envoyees dans la chimie: double precision
-
-      integer i,j,l,ic
-      REAL  temp_c(klev),press_c(klev)     ! T,p(mbar) a 1 lat donnee
-      REAL  declin_c                       ! declinaison en degres
-      REAL  cqy(klev,NC)                   ! frac mol qui seront modifiees
-      REAL  surfhaze(klev)
-      REAL  cprodaer(klev,4),cmaer(klev,4),ccsn(klev,4),ccsh(klev,4)
-      REAL  rinter(klev),nb(klev)
-      REAL  a,b
-      character str1*1,str2*2
-      integer ntotftop
-      parameter (ntotftop=50)
-      integer nftop(ntotftop),isaisonflux
-      REAL  fluxtop(NC),latit,tabletmp(ntotftop),factflux
-      character*10 nomsftop(ntotftop+1)
-      character*20 formt1,formt2,formt3
-      
-c     variables locales initialisees au premier appel
-
-      LOGICAL firstcal
-      DATA    firstcal/.true./
-      SAVE    firstcal
-      
-      REAL  mass(NC),duree
-      REAL  tablefluxtop(NC,jjm+1,5)
-      REAL  botCH4
-      DATA  botCH4/0.05/
-      REAL  krpd(15,ND+1,klev,jjm+1),krate(klev,NR)
-      integer reactif(5,NR),nom_prod(NC),nom_perte(NC)
-      integer prod(200,NC),perte(2,200,NC)
-      SAVE    mass,tablefluxtop,botCH4
-      SAVE    krpd,krate
-      SAVE    reactif,nom_prod,nom_perte,prod,perte
-      
-c-----------------------------------------------------------------------
-c***********************************************************************
-c
-c    Initialisations :
-c    ----------------
-
-      duree = dtchim  ! passage en real*4 pour appel a routines C
-
-      IF (firstcal) THEN
-            
-        print*,'CHIMIE, premier appel'
-        
-c ************************************
-c Au premier appel, initialisation de toutes les variables 
-c pour les routines de la chimie.
-c ************************************
-
-c Verification dimension verticale: coherence titan_for.h et klev
-c --------------------------------
-
-        if (klev.ne.NLEV) then
-	   print*,'PROBLEME de coherence dans la dimension verticale:'
-     .           ,klev,NLEV
-           STOP
-	endif
-
-c Verification du nombre de composes: coherence titan_for.h et nqmax-nmicro
-c ----------------------------------
-
-        if (ny.ne.NC) then
-	   print*,'PROBLEME de coherence dans le nombre de composes:'
-     .           ,ny,NC
-           STOP
-	endif
-
-c calcul de temp_c, densites et press_c a l'equateur:
-c --------------------------------------------------
-
-        print*,'pression, densites et temp a l equateur (chimie):'
-        print*,'level, press_c, nb, temp_c'
-        DO l=1,klev
-c     temp_c (K):
-         temp_c(l)  = ctemp(jjm/2+1,l)
-c     press_c (mbar):
-         press_c(l) = cplay(jjm/2+1,l)/100.
-c     nb (cm-3):
-         nb(l) = 1.e-4*press_c(l) / (RKBOL*temp_c(l))
-         print*,l,press_c(l),nb(l),temp_c(l)
-        ENDDO
-        
-c caracteristiques des composes:        
-c -----------------------------
-        mass = 0.0
-        call comp(nomqy_c,mass)
-        print*,'           Mass'
-        do ic=1,NC
-          print*,nomqy_c(ic),mass(ic)
-        enddo
-	
-
-c FLUX AU SOMMET: DEPENDANT DE LA LATITUDE ET DE LA SAISON... 
-
-c flux dans la derniere couche:(issu du modele 1D, a 500 km) 
-c -----------------------------
-c       !!  en cm-2.s-1 !!
-c  ces valeurs sont converties en cm-3.s-1 (dni/dt) dans a couche 54 
-c 
-c Lecture des tables de flux en fonction lat et saison
-c ----------------------------------------------------
-
-	WRITE(str2(1:2),'(i2.2)') ntotftop
-        formt1 = '(A10,'//str2//'(3X,A10))'
-        formt2 = '(F12.6,'//str2//'(2X,F13.6))'
-        formt3 = '(F13.6,'//str2//'(2X,F13.6))'
-	
-	do j=1,jjp1
-	  do ic=1,NC
-	    do l=1,5
-	     tablefluxtop(ic,j,l) = 0.
-	    enddo
-	  enddo
-	enddo
-	
-        do l=1,4
-	  WRITE(str1(1:1),'(i1.1)') l
-c hx1 -> Ls=RPI/2 / hx4 -> Ls=0
-c         open(10,file="../../INPUT/FLUXTOP/flux500.hs"//str1)
-          open(10,file="FLUXTOP/flux500.hs"//str1)
-c on remet l=1 et 5 -> Ls=0 / l=2 -> Ls=RPI/2 ...
-
-          read(10,formt1) nomsftop
-	  do j=1,ntotftop
-           do i=1,10   !justification a gauche...
-             if (nomsftop(j+1)(i:i).ne.' ') then
-                nomsftop(j) = nomsftop(j+1)(i:)
-                goto 100
-             endif
-           enddo
-100        continue
-c         print*,nomsftop(j)
-           nftop(j) = 0
-           do i=1,NC
-	    if (nomqy_c(i).eq.nomsftop(j)) then
-             nftop(j) = i
-	    endif
-	   enddo
-c          if (nftop(j).ne.0) print*,nomsftop(j),nomqy_c(nftop(j))
-	  enddo
-
-	 
-c lecture des flux. Les formats sont un peu alambiques... 
-c     a revoir eventuellement
-	  do j=1,jjp1/2+1
-	    read(10,formt2) latit,tabletmp
-	    do i=1,ntotftop
-              if (nftop(i).ne.0) then
-	       tablefluxtop(nftop(i),j,l+1) = tabletmp(i)
-              endif
-	    enddo
-	  enddo
-	  do j=jjp1/2+2,jjp1
-	    read(10,formt3) latit,tabletmp
-	    do i=1,ntotftop
-              if (nftop(i).ne.0) then
-	       tablefluxtop(nftop(i),j,l+1) = tabletmp(i)
-              endif
-	    enddo
-	  enddo
-
-	enddo  ! l
-	
-	do j=1,jjp1
-	  do ic=1,NC
-	     tablefluxtop(ic,j,1) = tablefluxtop(ic,j,5)
-	  enddo
-	enddo
-	
-c  Stockage des composes utilises dans la prod d'aerosols
-c     (aerprod=1) et dans H-> H2 (htoh2=1): utilaer
-c     ! decalage de 1 car utilise dans le c !
-c ------------------------------------------
-c + modifs du flux en fonction des obs UVIS et CIRS: C2H2,C2H6,HCN
-c + modifs des flux qui presentent un probleme evident: C3H8 et C4H10
-c ------------------------------------------
-
-        do ic=1,NC
-
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-c!!!remise de CH4 a 1.5%!!!!!!!!!!!!!!!!!!!!!!
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-c         if (nomqy_c(ic).eq."CH4") then
-c           do l=1,llm
-c            do j=1,ip1jmp1
-c              if (qy_c(j,l,ic).le.0.015) qy_c(j,l,ic) = 0.015 
-c            enddo
-c           enddo
-c         endif
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-          
-	  if (nomqy_c(ic).eq."C4H2") then
-            utilaer(10) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."HCN") then
-            utilaer(6)  = ic-1
-	    do j=1,jjp1
-              do i=1,5
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 6.
-              enddo
-            enddo
-	  endif
-	  if (nomqy_c(ic).eq."HC3N") then
-            utilaer(7)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."NCCN") then
-            utilaer(14) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."CH3CN") then
-            utilaer(15) = ic-1
-            utilaer(16) = ic-1 ! si pas C2H3CN, CH3CN utilise, mais reac. nulle
-	  endif
-	  if (nomqy_c(ic).eq."H") then
-            utilaer(1)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."H2") then
-            utilaer(2)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C2H2") then
-            utilaer(3)  = ic-1
-	    do j=1,jjp1
-              do i=1,5
-c             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 7.3
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 10.
-              enddo
-            enddo
-	  endif
-	  if (nomqy_c(ic).eq."AC6H6") then
-            utilaer(13) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C2H3CN") then
-            utilaer(16) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C2") then
-            utilaer(4)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C2H") then
-            utilaer(5)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C3N") then
-            utilaer(8)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."H2CN") then
-            utilaer(9)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C4H3") then
-            utilaer(11) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."AC6H5") then
-            utilaer(12) = ic-1
-	  endif
-
-	  if (nomqy_c(ic).eq."C2H6") then
-	    do j=1,jjp1
-              do i=1,5
-c             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 640.
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 1200.
-              enddo
-            enddo
-	  endif
-	  if (nomqy_c(ic).eq."C3H8") then
-	    do j=1,jjp1
-              do i=1,5
-c             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * (-1.)
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * (-3000.)
-              enddo
-            enddo
-	  endif
-	  if (nomqy_c(ic).eq."C4H10") then
-	    do j=1,jjp1
-              do i=1,5
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * (-1.)
-              enddo
-            enddo
-	  endif
-
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-c         if ((nomqy_c(ic).eq."HC3N").or.
-c    $        (nomqy_c(ic).eq."C3N")) then
-c     DO j=1,ip1jmp1 
-c       do l=1,34  ! p>~ 1 mbar
-c         qy_c(j,l,ic) = 1.e-30
-c       enddo
-c     ENDDO
-c         endif
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-        enddo
-	        
-c taux de photodissociations:
-c --------------------------
-        call disso(krpd,jjp1) 
-
-c reactions chimiques:
-c -------------------
-        call chimie(nomqy_c,nb,temp_c,krate,reactif,
-     .               nom_perte,nom_prod,perte,prod)
-c        print*,'nom_prod, nom_perte:'
-c        do ic=1,NC
-c          print*,nom_prod(ic),nom_perte(ic)
-c        enddo
-c        print*,'premieres prod, perte(1:reaction,2:compagnon):'
-c        do ic=1,NC
-c          print*,prod(1,ic),perte(1,1,ic),perte(2,1,ic)
-c        enddo
-
-c       l = klev-3
-c       print*,'krate a p=',press_c(l),' reactifs et produits:'
-c       do ic=1,ND+1
-c         print*,ic,krpd(7,ic,l,4)*nb(l),"  ",
-c    .     nomqy_c(reactif(1,ic)+1),
-c    .     nomqy_c(reactif(2,ic)+1),nomqy_c(reactif(3,ic)+1),
-c    .     nomqy_c(reactif(4,ic)+1),nomqy_c(reactif(5,ic)+1)
-c       enddo
-c       do ic=ND+2,NR
-c         print*,ic,krate(l,ic),"  ",
-c    .     nomqy_c(reactif(1,ic)+1),
-c    .     nomqy_c(reactif(2,ic)+1),nomqy_c(reactif(3,ic)+1),
-c    .     nomqy_c(reactif(4,ic)+1),nomqy_c(reactif(5,ic)+1)
-c       enddo
-
-      ENDIF  ! premier appel
-
-c***********************************************************************
-c-----------------------------------------------------------------------
-
-c   calcul declin_c (en degres)
-c   ---------------------------
-
-      declin_c = declin_rad*180./RPI
-c      print*,'declinaison en degre=',declin_c
-       
-c-----------------------------------------------------------------------
-c
-c   calcul du facteur d'interpolation entre deux saisons pour flux
-c   --------------------------------------------------------------
-
-       isaisonflux = int(ls_rad*2./RPI)+1
-       if (isaisonflux.eq.5) then
-         isaisonflux = 1
-       endif
-       factflux =    (sin(ls_rad)-sin((isaisonflux-1)*RPI/2.))/
-     .   (sin(isaisonflux*RPI/2.)-sin((isaisonflux-1)*RPI/2.))
-
-c***********************************************************************
-c***********************************************************************
-c
-c                BOUCLE SUR LES LATITUDES
-c
-      DO j=1,jjp1 
-      
-c***********************************************************************
-c***********************************************************************
-
-c-----------------------------------------------------------------------
-c
-c   Temperature, pression (mbar), densite (cm-3)
-c   -------------------------------------------
-
-       DO l=1,klev
-c     temp_c (K):
-               temp_c(l)  = ctemp(j,l)
-c     press_c (mbar):
-               press_c(l) = cplay(j,l)/100.
-c     nb (cm-3):
-               nb(l) = 1.e-4*press_c(l) / (RKBOL*temp_c(l))
-       ENDDO
-       
-c-----------------------------------------------------------------------
-c
-c   Distances radiales (intercouches, en km)
-c   ----------------------------------------
-
-       rinter(1) = RA/1000.  
-       do l=2,klev
-c A REVOIR: g doit varier avec r !
-         rinter(l) = rinter(l-1) + 
-     .    (cplev(j,l-1)-cplev(j,l))/cplay(j,l)*(RD*temp_c(l)/RG)/1000.
-c        print*,rinter(l)
-       enddo
-
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-c FLUX AU TOP: interpolation en saison et
-c conversion en cm-3 s-1 dans la couche klev-1 (NLEV-2 dans le C) 
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-        do ic=1,NC
-          fluxtop(ic) = tablefluxtop(ic,j,isaisonflux)  *(1-factflux)
-     .                + tablefluxtop(ic,j,isaisonflux+1)*  factflux
-          fluxtop(ic) = fluxtop(ic)/((rinter(klev)-rinter(klev-1))*1e5)
-        enddo
-
-c TEST: sortie de l'un des flux avec saveLs et factflux
-c       dans un fichier special pour tracer avec gnuplot
-c       if (j.eq.2) then
-c       open(unit=11,file="flux_80N.txt",status='old',position='append')
-c         write(11,*) ls_rad*180./RPI, factflux, 
-c    .           ((rinter(llm)-rinter(llm-1))*1e5), fluxtop
-c         write(11,*) ls_rad*180./RPI, factflux, 
-c    .     fluxtop(utilaer(3)+1)*((rinter(llm)-rinter(llm-1))*1e5),  !C2H2
-c    .     fluxtop(utilaer(6)+1)*((rinter(llm)-rinter(llm-1))*1e5)   !HCN
-c       close(11)
-c       endif
-c       if (j.eq.17) then
-c       open(unit=11,file="flux_eq.txt",status='old',position='append')
-c         write(11,*) ls_rad*180./RPI, factflux, 
-c    .           ((rinter(llm)-rinter(llm-1))*1e5), fluxtop
-c         write(11,*) ls_rad*180./RPI, factflux, 
-c    .     fluxtop(utilaer(3)+1)*((rinter(llm)-rinter(llm-1))*1e5),  !C2H2
-c    .     fluxtop(utilaer(6)+1)*((rinter(llm)-rinter(llm-1))*1e5)   !HCN
-c       close(11)
-c       endif
-c FIN TEST: sortie de l'un des flux avec ls et factflux
-
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-       if (firstcal.and.(j.eq.1)) then
-         print*,'Alt, densites et temp au pole (chimie):'
-         print*,'level, z_bas, nb, temp_c'
-         do l=1,klev
-          print*,l,rinter(l)-RA/1000.,nb(l),temp_c(l)
-         enddo
-       endif
-
-       if (firstcal.and.(j.eq.jjm/2)) then
-c         print*,'g,mugaz'
-c         print*,g,mugaz
-         print*,'Alt, densites et temp a l equateur (chimie):'
-         print*,'level, z_bas, nb, temp_c'
-         do l=1,klev
-          print*,l,rinter(l)-RA/1000.,nb(l),temp_c(l)
-         enddo
-       endif
-       
-c-----------------------------------------------------------------------
-c
-c   composition 
-c   ------------
-       
-       do ic=1,NC
-        do l=1,klev
-          cqy(l,ic) = qy_c(j,l,ic) 
-        enddo
-       enddo
-              
-c-----------------------------------------------------------------------
-c
-c   total haze area (um2/cm3)
-c   -------------------------
-
-       if (htoh2.eq.1) then
-        do l=1,klev
-! ATTENTION, INVERSION PAR RAPPORT A pg2.F !!!
-         surfhaze(l) = psurfhaze(j,klev+1-l) 
-c        if (j.eq.25)
-c    .    print*,'psurfhaze en um2/cm3:',surfhaze(l)
-        enddo
-       endif
-
-c-----------------------------------------------------------------------
-c
-c   Appel de chimietitan
-c   --------------------
-       
-       call gptitan(jjp1,rinter,temp_c,nb,
-     $              nomqy_c,cqy,fluxtop,
-     $              declin_c,duree,(j-1),mass,
-     $              botCH4,krpd,krate,reactif,
-     $              nom_prod,nom_perte,prod,perte,
-     $              aerprod,utilaer,cmaer,cprodaer,ccsn,ccsh,
-     $              htoh2,surfhaze)
-       
-c       if ( j.eq.jjm/2 ) 
-c    $        print*,cqy(1,1),cqy(klev,1),cqy(1,2),cqy(klev,2)
-c       if ( j.eq.jjm/2 ) 
-c    $  print*,qy_c(j,1,1),qy_c(j,klev,1),qy_c(j,1,2),qy_c(j,klev,2)
-
-c       stop
-
-c   Tendances composition 
-c   ---------------------
-     
-       do ic=1,NC
-        do l=1,klev
-          dqyc(j,l,ic) = (cqy(l,ic) - qy_c(j,l,ic))/dtchim  ! en /s
-        enddo
-       enddo
-
-c-----------------------------------------------------------------------
-c
-c   production aer 
-c   --------------
-       
-       if (aerprod.eq.1) then
-
-       do ic=1,4
-        do l=1,klev
-          prodaer(j,l,ic) = cprodaer(l,ic) 
-          maer(j,l,ic)    = cmaer(l,ic) 
-          csn(j,l,ic)     = ccsn(l,ic) 
-          csh(j,l,ic)     = ccsh(l,ic) 
-        enddo
-       enddo
-
-       endif
-
-c***********************************************************************
-c***********************************************************************
-c
-c              FIN: BOUCLE SUR LES LATITUDES 
-c
-      ENDDO
-     
-c***********************************************************************
-c***********************************************************************
-
-
-      firstcal = .false.
-      RETURN
-      END
Index: unk/libf/phytitan/calchim.spleprec
===================================================================
--- /trunk/libf/phytitan/calchim.spleprec	(revision 101)
+++ 	(revision )
@@ -1,565 +1,0 @@
-      SUBROUTINE calchim(ny,qy_c,nomqy_c,declin_rad,ls_rad,dtchim,
-     .                   ctemp,cplay,cplev,
-     .                   
-     .
-     .                   dqyc)
-      implicit none
-      
-c-------------------------------------------------
-c  
-c     Introduction d'une routine chimique
-c
-c     Auteur: S. Lebonnois,  01/2000 | 09/2003 
-c            adaptation pour Titan 3D: 02/2009
-c
-c-------------------------------------------------
-c
-
-#include "dimensions.h"
-#include "dimphy.h"
-#include "clesphys.h"
-#include "paramet.h"
-#include "YOMCST.h"
-
-#include "titan_for.h"  
-!!!  doit etre en accord avec titan.h
-#include "aerprod.h"
-
-c    Arguments
-c    ---------
-
-      INTEGER      ny                     ! nb de composes (nqmax-nmicro) 
-      REAL         qy_c(jjm+1,klev,NC)    ! Especes chimiques apres adv.+diss.
-      character*10 nomqy_c(NC+1)          ! Noms des especes chimiques
-      REAL         declin_rad,ls_rad      ! declinaison et long solaire en radians
-      REAL         dtchim                 ! pas de temps chimie
-      REAL         ctemp(jjm+1,klev)      ! Temperature
-      REAL         cplay(jjm+1,klev)      ! pression (Pa)
-      REAL         cplev(jjm+1,klev)      ! pression intercouches (Pa)
-      
-      REAL         dqyc(jjm+1,klev,nqmx)  ! Tendances especes chimiques (nqmx, mais en fait NC...)
-      
-c    Local variables :
-c    -----------------
-c variables envoyees dans la chimie: simple precision
-
-      integer i,j,l,ic
-      REAL*4  temp_c(klev),press_c(klev)     ! T,p(mbar) a 1 lat donnee
-      REAL*4  declin_c                       ! declinaison en degres
-      REAL*4  cqy(klev,NC)                   ! frac mol qui seront modifiees
-      REAL*4  surfhaze(klev)
-      REAL*4  cprodaer(klev,4),cmaer(klev,4),ccsn(klev,4),ccsh(klev,4)
-      REAL*4  rinter(klev),nb(klev)
-      REAL*4  a,b
-      character str1*1,str2*2
-      integer ntotftop
-      parameter (ntotftop=50)
-      integer nftop(ntotftop),isaisonflux
-      REAL*4  fluxtop(NC),latit,tabletmp(ntotftop),factflux
-      character*10 nomsftop(ntotftop+1)
-      character*20 formt1,formt2,formt3
-      
-c     variables locales initialisees au premier appel
-
-      LOGICAL firstcal
-      DATA    firstcal/.true./
-      SAVE    firstcal
-      
-      REAL*4  mass(NC),duree
-      REAL*4  tablefluxtop(NC,jjm+1,5)
-      REAL*4  botCH4
-      DATA    botCH4/0.05/
-      REAL*4  krpd(15,ND+1,klev,jjm+1),krate(klev,NR)
-      integer reactif(5,NR),nom_prod(NC),nom_perte(NC)
-      integer prod(200,NC),perte(2,200,NC)
-      SAVE    mass,tablefluxtop,botCH4
-      SAVE    krpd,krate
-      SAVE    reactif,nom_prod,nom_perte,prod,perte
-      
-c-----------------------------------------------------------------------
-c***********************************************************************
-c
-c    Initialisations :
-c    ----------------
-
-      duree = dtchim  ! passage en real*4 pour appel a routines C
-
-      IF (firstcal) THEN
-            
-        print*,'CHIMIE, premier appel'
-        
-c ************************************
-c Au premier appel, initialisation de toutes les variables 
-c pour les routines de la chimie.
-c ************************************
-
-c Verification dimension verticale: coherence titan_for.h et klev
-c --------------------------------
-
-        if (klev.ne.NLEV) then
-	   print*,'PROBLEME de coherence dans la dimension verticale:'
-     .           ,klev,NLEV
-           STOP
-	endif
-
-c Verification du nombre de composes: coherence titan_for.h et nqmax-nmicro
-c ----------------------------------
-
-        if (ny.ne.NC) then
-	   print*,'PROBLEME de coherence dans le nombre de composes:'
-     .           ,ny,NC
-           STOP
-	endif
-
-c calcul de temp_c, densites et press_c a l'equateur:
-c --------------------------------------------------
-
-        print*,'pression, densites et temp a l equateur (chimie):'
-        print*,'level, press_c, nb, temp_c'
-        DO l=1,klev
-c     temp_c (K):
-         temp_c(l)  = ctemp(jjm/2+1,l)
-c     press_c (mbar):
-         press_c(l) = cplay(jjm/2+1,l)/100.
-c     nb (cm-3):
-         nb(l) = 1.e-4*press_c(l) / (RKBOL*temp_c(l))
-         print*,l,press_c(l),nb(l),temp_c(l)
-        ENDDO
-        
-c caracteristiques des composes:        
-c -----------------------------
-        mass = 0.0
-        call comp(nomqy_c,mass)
-        print*,'           Mass'
-        do ic=1,NC
-          print*,nomqy_c(ic),mass(ic)
-        enddo
-	
-
-c FLUX AU SOMMET: DEPENDANT DE LA LATITUDE ET DE LA SAISON... 
-
-c flux dans la derniere couche:(issu du modele 1D, a 500 km) 
-c -----------------------------
-c       !!  en cm-2.s-1 !!
-c  ces valeurs sont converties en cm-3.s-1 (dni/dt) dans a couche 54 
-c 
-c Lecture des tables de flux en fonction lat et saison
-c ----------------------------------------------------
-
-	WRITE(str2(1:2),'(i2.2)') ntotftop
-        formt1 = '(A10,'//str2//'(3X,A10))'
-        formt2 = '(F12.6,'//str2//'(2X,F13.6))'
-        formt3 = '(F13.6,'//str2//'(2X,F13.6))'
-	
-	do j=1,jjp1
-	  do ic=1,NC
-	    do l=1,5
-	     tablefluxtop(ic,j,l) = 0.
-	    enddo
-	  enddo
-	enddo
-	
-        do l=1,4
-	  WRITE(str1(1:1),'(i1.1)') l
-c hx1 -> Ls=RPI/2 / hx4 -> Ls=0
-c         open(10,file="../../INPUT/FLUXTOP/flux500.hs"//str1)
-          open(10,file="FLUXTOP/flux500.hs"//str1)
-c on remet l=1 et 5 -> Ls=0 / l=2 -> Ls=RPI/2 ...
-
-          read(10,formt1) nomsftop
-	  do j=1,ntotftop
-           do i=1,10   !justification a gauche...
-             if (nomsftop(j+1)(i:i).ne.' ') then
-                nomsftop(j) = nomsftop(j+1)(i:)
-                goto 100
-             endif
-           enddo
-100        continue
-c         print*,nomsftop(j)
-           nftop(j) = 0
-           do i=1,NC
-	    if (nomqy_c(i).eq.nomsftop(j)) then
-             nftop(j) = i
-	    endif
-	   enddo
-c          if (nftop(j).ne.0) print*,nomsftop(j),nomqy_c(nftop(j))
-	  enddo
-
-	 
-c lecture des flux. Les formats sont un peu alambiques... 
-c     a revoir eventuellement
-	  do j=1,jjp1/2+1
-	    read(10,formt2) latit,tabletmp
-	    do i=1,ntotftop
-              if (nftop(i).ne.0) then
-	       tablefluxtop(nftop(i),j,l+1) = tabletmp(i)
-              endif
-	    enddo
-	  enddo
-	  do j=jjp1/2+2,jjp1
-	    read(10,formt3) latit,tabletmp
-	    do i=1,ntotftop
-              if (nftop(i).ne.0) then
-	       tablefluxtop(nftop(i),j,l+1) = tabletmp(i)
-              endif
-	    enddo
-	  enddo
-
-	enddo  ! l
-	
-	do j=1,jjp1
-	  do ic=1,NC
-	     tablefluxtop(ic,j,1) = tablefluxtop(ic,j,5)
-	  enddo
-	enddo
-	
-c  Stockage des composes utilises dans la prod d'aerosols
-c     (aerprod=1) et dans H-> H2 (htoh2=1): utilaer
-c     ! decalage de 1 car utilise dans le c !
-c ------------------------------------------
-c + modifs du flux en fonction des obs UVIS et CIRS: C2H2,C2H6,HCN
-c + modifs des flux qui presentent un probleme evident: C3H8 et C4H10
-c ------------------------------------------
-
-        do ic=1,NC
-
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-c!!!remise de CH4 a 1.5%!!!!!!!!!!!!!!!!!!!!!!
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-c         if (nomqy_c(ic).eq."CH4") then
-c           do l=1,llm
-c            do j=1,ip1jmp1
-c              if (qy_c(j,l,ic).le.0.015) qy_c(j,l,ic) = 0.015 
-c            enddo
-c           enddo
-c         endif
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-          
-	  if (nomqy_c(ic).eq."C4H2") then
-            utilaer(10) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."HCN") then
-            utilaer(6)  = ic-1
-	    do j=1,jjp1
-              do i=1,5
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 6.
-              enddo
-            enddo
-	  endif
-	  if (nomqy_c(ic).eq."HC3N") then
-            utilaer(7)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."NCCN") then
-            utilaer(14) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."CH3CN") then
-            utilaer(15) = ic-1
-            utilaer(16) = ic-1 ! si pas C2H3CN, CH3CN utilise, mais reac. nulle
-	  endif
-	  if (nomqy_c(ic).eq."H") then
-            utilaer(1)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."H2") then
-            utilaer(2)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C2H2") then
-            utilaer(3)  = ic-1
-	    do j=1,jjp1
-              do i=1,5
-c             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 7.3
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 10.
-              enddo
-            enddo
-	  endif
-	  if (nomqy_c(ic).eq."AC6H6") then
-            utilaer(13) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C2H3CN") then
-            utilaer(16) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C2") then
-            utilaer(4)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C2H") then
-            utilaer(5)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C3N") then
-            utilaer(8)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."H2CN") then
-            utilaer(9)  = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."C4H3") then
-            utilaer(11) = ic-1
-	  endif
-	  if (nomqy_c(ic).eq."AC6H5") then
-            utilaer(12) = ic-1
-	  endif
-
-	  if (nomqy_c(ic).eq."C2H6") then
-	    do j=1,jjp1
-              do i=1,5
-c             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 640.
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * 1200.
-              enddo
-            enddo
-	  endif
-	  if (nomqy_c(ic).eq."C3H8") then
-	    do j=1,jjp1
-              do i=1,5
-c             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * (-1.)
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * (-3000.)
-              enddo
-            enddo
-	  endif
-	  if (nomqy_c(ic).eq."C4H10") then
-	    do j=1,jjp1
-              do i=1,5
-             tablefluxtop(ic,j,i) = tablefluxtop(ic,j,i) * (-1.)
-              enddo
-            enddo
-	  endif
-
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-c         if ((nomqy_c(ic).eq."HC3N").or.
-c    $        (nomqy_c(ic).eq."C3N")) then
-c     DO j=1,ip1jmp1 
-c       do l=1,34  ! p>~ 1 mbar
-c         qy_c(j,l,ic) = 1.e-30
-c       enddo
-c     ENDDO
-c         endif
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-        enddo
-	        
-c taux de photodissociations:
-c --------------------------
-        call disso(krpd,jjp1) 
-
-c reactions chimiques:
-c -------------------
-        call chimie(nomqy_c,nb,temp_c,krate,reactif,
-     .               nom_perte,nom_prod,perte,prod)
-c        print*,'nom_prod, nom_perte:'
-c        do ic=1,NC
-c          print*,nom_prod(ic),nom_perte(ic)
-c        enddo
-c        print*,'premieres prod, perte(1:reaction,2:compagnon):'
-c        do ic=1,NC
-c          print*,prod(1,ic),perte(1,1,ic),perte(2,1,ic)
-c        enddo
-
-c       l = klev-3
-c       print*,'krate a p=',press_c(l),' reactifs et produits:'
-c       do ic=1,ND+1
-c         print*,ic,krpd(7,ic,l,4)*nb(l),"  ",
-c    .     nomqy_c(reactif(1,ic)+1),
-c    .     nomqy_c(reactif(2,ic)+1),nomqy_c(reactif(3,ic)+1),
-c    .     nomqy_c(reactif(4,ic)+1),nomqy_c(reactif(5,ic)+1)
-c       enddo
-c       do ic=ND+2,NR
-c         print*,ic,krate(l,ic),"  ",
-c    .     nomqy_c(reactif(1,ic)+1),
-c    .     nomqy_c(reactif(2,ic)+1),nomqy_c(reactif(3,ic)+1),
-c    .     nomqy_c(reactif(4,ic)+1),nomqy_c(reactif(5,ic)+1)
-c       enddo
-
-      ENDIF  ! premier appel
-
-c***********************************************************************
-c-----------------------------------------------------------------------
-
-c   calcul declin_c (en degres)
-c   ---------------------------
-
-      declin_c = declin_rad*180./RPI
-c      print*,'declinaison en degre=',declin_c
-       
-c-----------------------------------------------------------------------
-c
-c   calcul du facteur d'interpolation entre deux saisons pour flux
-c   --------------------------------------------------------------
-
-       isaisonflux = int(ls_rad*2./RPI)+1
-       if (isaisonflux.eq.5) then
-         isaisonflux = 1
-       endif
-       factflux =    (sin(ls_rad)-sin((isaisonflux-1)*RPI/2.))/
-     .   (sin(isaisonflux*RPI/2.)-sin((isaisonflux-1)*RPI/2.))
-
-c***********************************************************************
-c***********************************************************************
-c
-c                BOUCLE SUR LES LATITUDES
-c
-      DO j=1,jjp1 
-      
-c***********************************************************************
-c***********************************************************************
-
-c-----------------------------------------------------------------------
-c
-c   Temperature, pression (mbar), densite (cm-3)
-c   -------------------------------------------
-
-       DO l=1,klev
-c     temp_c (K):
-               temp_c(l)  = ctemp(j,l)
-c     press_c (mbar):
-               press_c(l) = cplay(j,l)/100.
-c     nb (cm-3):
-               nb(l) = 1.e-4*press_c(l) / (RKBOL*temp_c(l))
-       ENDDO
-       
-c-----------------------------------------------------------------------
-c
-c   Distances radiales (intercouches, en km)
-c   ----------------------------------------
-
-       rinter(1) = RA/1000.  
-       do l=2,klev
-c A REVOIR: g doit varier avec r !
-         rinter(l) = rinter(l-1) + 
-     .    (cplev(j,l-1)-cplev(j,l))/cplay(j,l)*(RD*temp_c(l)/RG)/1000.
-c        print*,rinter(l)
-       enddo
-
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-c FLUX AU TOP: interpolation en saison et
-c conversion en cm-3 s-1 dans la couche klev-1 (NLEV-2 dans le C) 
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-        do ic=1,NC
-          fluxtop(ic) = tablefluxtop(ic,j,isaisonflux)  *(1-factflux)
-     .                + tablefluxtop(ic,j,isaisonflux+1)*  factflux
-          fluxtop(ic) = fluxtop(ic)/((rinter(klev)-rinter(klev-1))*1e5)
-        enddo
-
-c TEST: sortie de l'un des flux avec saveLs et factflux
-c       dans un fichier special pour tracer avec gnuplot
-c       if (j.eq.2) then
-c       open(unit=11,file="flux_80N.txt",status='old',position='append')
-c         write(11,*) ls_rad*180./RPI, factflux, 
-c    .           ((rinter(llm)-rinter(llm-1))*1e5), fluxtop
-c         write(11,*) ls_rad*180./RPI, factflux, 
-c    .     fluxtop(utilaer(3)+1)*((rinter(llm)-rinter(llm-1))*1e5),  !C2H2
-c    .     fluxtop(utilaer(6)+1)*((rinter(llm)-rinter(llm-1))*1e5)   !HCN
-c       close(11)
-c       endif
-c       if (j.eq.17) then
-c       open(unit=11,file="flux_eq.txt",status='old',position='append')
-c         write(11,*) ls_rad*180./RPI, factflux, 
-c    .           ((rinter(llm)-rinter(llm-1))*1e5), fluxtop
-c         write(11,*) ls_rad*180./RPI, factflux, 
-c    .     fluxtop(utilaer(3)+1)*((rinter(llm)-rinter(llm-1))*1e5),  !C2H2
-c    .     fluxtop(utilaer(6)+1)*((rinter(llm)-rinter(llm-1))*1e5)   !HCN
-c       close(11)
-c       endif
-c FIN TEST: sortie de l'un des flux avec ls et factflux
-
-c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-       if (firstcal.and.(j.eq.1)) then
-         print*,'Alt, densites et temp au pole (chimie):'
-         print*,'level, z_bas, nb, temp_c'
-         do l=1,klev
-          print*,l,rinter(l)-RA/1000.,nb(l),temp_c(l)
-         enddo
-       endif
-
-       if (firstcal.and.(j.eq.jjm/2)) then
-c         print*,'g,mugaz'
-c         print*,g,mugaz
-         print*,'Alt, densites et temp a l equateur (chimie):'
-         print*,'level, z_bas, nb, temp_c'
-         do l=1,klev
-          print*,l,rinter(l)-RA/1000.,nb(l),temp_c(l)
-         enddo
-       endif
-       
-c-----------------------------------------------------------------------
-c
-c   composition 
-c   ------------
-       
-       do ic=1,NC
-        do l=1,klev
-          cqy(l,ic) = qy_c(j,l,ic) 
-        enddo
-       enddo
-              
-c-----------------------------------------------------------------------
-c
-c   total haze area (um2/cm3)
-c   -------------------------
-
-       if (htoh2.eq.1) then
-        do l=1,klev
-! ATTENTION, INVERSION PAR RAPPORT A pg2.F !!!
-         surfhaze(l) = psurfhaze(j,klev+1-l) 
-c        if (j.eq.25)
-c    .    print*,'psurfhaze en um2/cm3:',surfhaze(l)
-        enddo
-       endif
-
-c-----------------------------------------------------------------------
-c
-c   Appel de chimietitan
-c   --------------------
-       
-       call gptitan(jjp1,rinter,temp_c,nb,
-     $              nomqy_c,cqy,fluxtop,
-     $              declin_c,duree,(j-1),mass,
-     $              botCH4,krpd,krate,reactif,
-     $              nom_prod,nom_perte,prod,perte,
-     $              aerprod,utilaer,cmaer,cprodaer,ccsn,ccsh,
-     $              htoh2,surfhaze)
-       
-c       if ( j.eq.jjm/2 ) 
-c    $        print*,cqy(1,1),cqy(klev,1),cqy(1,2),cqy(klev,2)
-c       if ( j.eq.jjm/2 ) 
-c    $  print*,qy_c(j,1,1),qy_c(j,klev,1),qy_c(j,1,2),qy_c(j,klev,2)
-
-c       stop
-
-c   Tendances composition 
-c   ---------------------
-     
-       do ic=1,NC
-        do l=1,klev
-          dqyc(j,l,ic) = (cqy(l,ic) - qy_c(j,l,ic))/dtchim  ! en /s
-        enddo
-       enddo
-
-c-----------------------------------------------------------------------
-c
-c   production aer 
-c   --------------
-       
-       if (aerprod.eq.1) then
-
-       do ic=1,4
-        do l=1,klev
-          prodaer(j,l,ic) = cprodaer(l,ic) 
-          maer(j,l,ic)    = cmaer(l,ic) 
-          csn(j,l,ic)     = ccsn(l,ic) 
-          csh(j,l,ic)     = ccsh(l,ic) 
-        enddo
-       enddo
-
-       endif
-
-c***********************************************************************
-c***********************************************************************
-c
-c              FIN: BOUCLE SUR LES LATITUDES 
-c
-      ENDDO
-     
-c***********************************************************************
-c***********************************************************************
-
-
-      firstcal = .false.
-      RETURN
-      END
Index: /trunk/libf/phytitan/clmain.F
===================================================================
--- /trunk/libf/phytitan/clmain.F	(revision 101)
+++ /trunk/libf/phytitan/clmain.F	(revision 102)
@@ -34,4 +34,5 @@
       USE ioipsl
       USE interface_surf
+      use dimphy
       IMPLICIT none
 c======================================================================
@@ -65,5 +66,4 @@
 c======================================================================
 #include "dimensions.h"
-#include "dimphy.h"
 c$$$ PB ajout pour soil
 #include "dimsoil.h"
@@ -477,4 +477,5 @@
 
       USE interface_surf
+      use dimphy
 
       IMPLICIT none
@@ -484,5 +485,4 @@
 c======================================================================
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 #include "dimsoil.h"
@@ -755,4 +755,6 @@
      e                  paprs,pplay,delp,
      s                  d_ven,flux_v)
+
+      use dimphy
       IMPLICIT none
 c======================================================================
@@ -778,5 +780,4 @@
 c======================================================================
 #include "dimensions.h"
-#include "dimphy.h"
 #include "iniprint.h"
       INTEGER knon
@@ -883,4 +884,6 @@
      .                  ts,u,v,t,
      .                  pcfm, pcfh)
+
+      use dimphy
       IMPLICIT none
 c======================================================================
@@ -903,5 +906,4 @@
 c======================================================================
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 #include "iniprint.h"
@@ -938,7 +940,4 @@
       PARAMETER (prandtl=0.4)
       INTEGER isommet ! le sommet de la couche limite
-c TEST VENUS
-      PARAMETER (isommet=klev)
-c     PARAMETER (isommet=5)
 
       LOGICAL tvirtu ! calculer Ri d'une maniere plus performante
@@ -984,4 +983,6 @@
       DATA appel1er /.TRUE./
 c
+      isommet=klev
+
       IF (appel1er) THEN
         if (prt_level > 9) THEN
@@ -1166,4 +1167,6 @@
       SUBROUTINE coefkz2(knon, paprs, pplay,t,
      .                  pcfm, pcfh)
+
+      use dimphy
       IMPLICIT none
 c======================================================================
@@ -1182,5 +1185,4 @@
 c======================================================================
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 #include "iniprint.h"
Index: /trunk/libf/phytitan/cltrac.F
===================================================================
--- /trunk/libf/phytitan/cltrac.F	(revision 101)
+++ /trunk/libf/phytitan/cltrac.F	(revision 102)
@@ -4,4 +4,6 @@
       SUBROUTINE cltrac(dtime,coef,t,tr,flux,paprs,pplay,delp,
      s                  d_tr)
+
+      use dimphy
       IMPLICIT none
 c======================================================================
@@ -26,5 +28,4 @@
 c======================================================================
 #include "dimensions.h"
-#include "dimphy.h"
       REAL dtime
       REAL coef(klon,klev)
@@ -118,4 +119,11 @@
       ENDDO
 c
+c ATTENTION SHUNTE!!!!!!
+c     DO k = 1, klev
+c     DO i = 1, klon
+c        d_tr(i,k) = 0.
+c     ENDDO
+c     ENDDO
+
       RETURN
       END
Index: /trunk/libf/phytitan/comgeomphy.F90
===================================================================
--- /trunk/libf/phytitan/comgeomphy.F90	(revision 102)
+++ /trunk/libf/phytitan/comgeomphy.F90	(revision 102)
@@ -0,0 +1,23 @@
+module comgeomphy
+   real,save,allocatable :: airephy(:)
+   real,save,allocatable :: cuphy(:)
+   real,save,allocatable :: cvphy(:)
+   real,save,allocatable :: rlatd(:)
+   real,save,allocatable :: rlond(:)
+!$OMP THREADPRIVATE(airephy,cuphy,cvphy,rlatd,rlond)
+contains
+  
+  subroutine initcomgeomphy
+  use dimphy
+  implicit none
+    
+ 
+    allocate(airephy(klon))
+    allocate(cuphy(klon))
+    allocate(cvphy(klon))
+    allocate(rlatd(klon))
+    allocate(rlond(klon))
+
+  end subroutine initcomgeomphy
+  
+end module comgeomphy
Index: unk/libf/phytitan/comgeomphy.h
===================================================================
--- /trunk/libf/phytitan/comgeomphy.h	(revision 101)
+++ 	(revision )
@@ -1,9 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/comgeomphy.h,v 1.1.1.1 2004/05/19 12:53:09 lmdzadmin Exp $
-!
-c
-c
-c Common de passage de la geometrie de la dynamique a la physique
-      real airephy(klon),cuphy(klon),cvphy(klon)
-      REAL rlatd(klon), rlond(klon)
-      common/comgeomphy/airephy,cuphy,cvphy,rlatd, rlond
Index: /trunk/libf/phytitan/cooling.F
===================================================================
--- /trunk/libf/phytitan/cooling.F	(revision 101)
+++ /trunk/libf/phytitan/cooling.F	(revision 102)
@@ -1,4 +1,3 @@
       SUBROUTINE COOLING(ngrid,NL,PRESS,TEMP,Z,Q0,lwnet,pfluxi)
-      IMPLICIT NONE
 
 c=======================================================================
@@ -52,6 +51,7 @@
 c   ------------
 
+      use dimphy
+      IMPLICIT NONE
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
       INTEGER NLAYER,NSPECI,NSPC1I
Index: /trunk/libf/phytitan/diagphy.F
===================================================================
--- /trunk/libf/phytitan/diagphy.F	(revision 101)
+++ /trunk/libf/phytitan/diagphy.F	(revision 102)
@@ -7,4 +7,7 @@
      $    , d_etp_tot, d_qt_tot, d_ec_tot
      $    , fs_bound, fq_bound)
+
+! ATTENTION !! PAS DU TOUT A JOUR POUR VENUS OU TITAN...
+
 C======================================================================
 C
@@ -45,8 +48,8 @@
 C======================================================================
 C 
+      use dimphy
       implicit none
 
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 C
@@ -205,8 +208,8 @@
 c======================================================================
  
+      use dimphy
       IMPLICIT NONE
 C
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 C
Index: /trunk/libf/phytitan/dimphy.F90
===================================================================
--- /trunk/libf/phytitan/dimphy.F90	(revision 102)
+++ /trunk/libf/phytitan/dimphy.F90	(revision 102)
@@ -0,0 +1,40 @@
+MODULE dimphy
+  
+  INTEGER,SAVE :: klon
+  INTEGER,SAVE :: kdlon
+  INTEGER,SAVE :: kfdia
+  INTEGER,SAVE :: kidia
+  INTEGER,SAVE :: klev
+  INTEGER,SAVE :: klevp1
+  INTEGER,SAVE :: klevm1
+  INTEGER,SAVE :: kflev
+
+!$OMP THREADPRIVATE(klon,kfdia,kidia,kdlon)
+  REAL,save,allocatable,dimension(:) :: zmasq
+!$OMP THREADPRIVATE(zmasq)   
+
+CONTAINS
+  
+  SUBROUTINE init_dimphy(klon0,klev0)
+  IMPLICIT NONE
+  
+    INTEGER, INTENT(in) :: klon0
+    INTEGER, INTENT(in) :: klev0
+    
+    klon=klon0
+    
+    kdlon=klon
+    kidia=1
+    kfdia=klon
+!$OMP MASTER 
+    klev=klev0
+    klevp1=klev+1
+    klevm1=klev-1
+    kflev=klev
+!$OMP END MASTER    
+    ALLOCATE(zmasq(klon))    
+    
+  END SUBROUTINE init_dimphy
+
+  
+END MODULE dimphy
Index: unk/libf/phytitan/dimphy.h
===================================================================
--- /trunk/libf/phytitan/dimphy.h	(revision 101)
+++ 	(revision )
@@ -1,7 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/dimphy.h,v 1.1.1.1 2004/05/19 12:53:07 lmdzadmin Exp $
-!
-c-----------------------------------------------------------------------
-      INTEGER KIDIA, KFDIA, KLON, KLEV
-      PARAMETER (KIDIA=1,KFDIA=iim*(jjm-1)+2-1/jjm,
-     .           KLON=KFDIA-KIDIA+1,KLEV=llm)
Index: /trunk/libf/phytitan/heating.F
===================================================================
--- /trunk/libf/phytitan/heating.F	(revision 101)
+++ /trunk/libf/phytitan/heating.F	(revision 102)
@@ -1,5 +1,3 @@
        SUBROUTINE heating(dist,rmu0,fract,sol_htg,swnet)
-
-      IMPLICIT NONE
 
 
@@ -30,6 +28,7 @@
 c   -------------
 
+      use dimphy
+      IMPLICIT NONE
 #include "dimensions.h"
-#include "dimphy.h"
 
       INTEGER NLEVEL,NLAYER,NSPECV
Index: /trunk/libf/phytitan/hgardfou.F
===================================================================
--- /trunk/libf/phytitan/hgardfou.F	(revision 101)
+++ /trunk/libf/phytitan/hgardfou.F	(revision 102)
@@ -3,10 +3,10 @@
 !
       SUBROUTINE hgardfou (t,tsol,text)
-      IMPLICIT none
 c======================================================================
 c Verifier la temperature
 c======================================================================
+      use dimphy
+      IMPLICIT none
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
       REAL t(klon,klev), tsol(klon)
Index: /trunk/libf/phytitan/inicondens.F
===================================================================
--- /trunk/libf/phytitan/inicondens.F	(revision 101)
+++ /trunk/libf/phytitan/inicondens.F	(revision 102)
@@ -9,7 +9,7 @@
 c   -------------
 
+      use dimphy
       IMPLICIT NONE
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 
Index: /trunk/libf/phytitan/iniphysiq.F
===================================================================
--- /trunk/libf/phytitan/iniphysiq.F	(revision 101)
+++ /trunk/libf/phytitan/iniphysiq.F	(revision 102)
@@ -9,5 +9,5 @@
      $           plat,plon,parea,pcu,pcv,
      $           prad,pg,pr,pcpp)
-      IMPLICIT NONE
+
 c
 c=======================================================================
@@ -43,7 +43,8 @@
 c   -------------
  
+      use dimphy
+      USE comgeomphy
+      IMPLICIT NONE
 #include "dimensions.h"
-#include "dimphy.h"
-#include "comgeomphy.h"
 
       REAL prad,pg,pr,pcpp,punjours
Index: /trunk/libf/phytitan/init_phys_lmdz.F90
===================================================================
--- /trunk/libf/phytitan/init_phys_lmdz.F90	(revision 102)
+++ /trunk/libf/phytitan/init_phys_lmdz.F90	(revision 102)
@@ -0,0 +1,19 @@
+!
+!$Header$
+!
+SUBROUTINE init_phys_lmdz(iim,jjp1,llm,nb_proc,distrib)
+  USE mod_grid_phy_lmdz
+  USE dimphy
+  IMPLICIT NONE
+  
+    INTEGER,INTENT(in) :: iim
+    INTEGER,INTENT(in) :: jjp1
+    INTEGER,INTENT(in) :: llm
+    INTEGER,INTENT(in) :: nb_proc
+    INTEGER,INTENT(in) :: distrib(0:nb_proc-1)
+
+
+    CALL init_grid_phy_lmdz(iim,jjp1,llm)
+    CALL init_dimphy(klon,nbp_lev)
+ 
+END SUBROUTINE init_phys_lmdz  
Index: /trunk/libf/phytitan/mod_grid_phy_lmdz.F90
===================================================================
--- /trunk/libf/phytitan/mod_grid_phy_lmdz.F90	(revision 102)
+++ /trunk/libf/phytitan/mod_grid_phy_lmdz.F90	(revision 102)
@@ -0,0 +1,447 @@
+!
+!$Header$
+!
+MODULE mod_grid_phy_lmdz
+  INTEGER,SAVE :: nbp_lon  ! == iim
+  INTEGER,SAVE :: nbp_lat  ! == jjmp1
+  INTEGER,SAVE :: nbp_lev  ! == llm
+  INTEGER,SAVE :: klon_glo
+
+  INTERFACE grid1dTo2d_glo
+    MODULE PROCEDURE grid1dTo2d_glo_i,grid1dTo2d_glo_i1,grid1dTo2d_glo_i2,grid1dTo2d_glo_i3, &
+                     grid1dTo2d_glo_r,grid1dTo2d_glo_r1,grid1dTo2d_glo_r2,grid1dTo2d_glo_r3, &
+		     grid1dTo2d_glo_l,grid1dTo2d_glo_l1,grid1dTo2d_glo_l2,grid1dTo2d_glo_l3
+   END INTERFACE 
+
+   INTERFACE grid2dTo1d_glo
+    MODULE PROCEDURE grid2dTo1d_glo_i,grid2dTo1d_glo_i1,grid2dTo1d_glo_i2,grid2dTo1d_glo_i3, &
+                     grid2dTo1d_glo_r,grid2dTo1d_glo_r1,grid2dTo1d_glo_r2,grid2dTo1d_glo_r3, &
+		     grid2dTo1d_glo_l,grid2dTo1d_glo_l1,grid2dTo1d_glo_l2,grid2dTo1d_glo_l3
+   END INTERFACE 
+ 
+CONTAINS
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+!! SUBROUTINE grid1dTo2d  !!  
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+
+  SUBROUTINE init_grid_phy_lmdz(iim,jjp1,llm)
+  IMPLICIT NONE
+  INTEGER, INTENT(in) :: iim
+  INTEGER, INTENT(in) :: jjp1
+  INTEGER, INTENT(in) :: llm
+  
+    nbp_lon=iim
+    nbp_lat=jjp1
+    nbp_lev=llm
+    klon_glo=(iim*jjp1)-2*(iim-1)
+  
+  END SUBROUTINE init_grid_phy_lmdz
+  
+  
+  SUBROUTINE grid1dTo2d_glo_i(VarIn,VarOut)  
+  IMPLICIT NONE  
+    INTEGER,INTENT(IN),DIMENSION(:)     :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_igen(VarIn,VarOut,1)
+  
+  END SUBROUTINE grid1dTo2d_glo_i
+  
+
+  SUBROUTINE grid1dTo2d_glo_i1(VarIn,VarOut)  
+  IMPLICIT NONE  
+    INTEGER,INTENT(IN),DIMENSION(:,:)     :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(:,:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_igen(VarIn,VarOut,size(VarIn,2))
+  
+  END SUBROUTINE grid1dTo2d_glo_i1
+
+  SUBROUTINE grid1dTo2d_glo_i2(VarIn,VarOut)  
+  IMPLICIT NONE  
+    INTEGER,INTENT(IN),DIMENSION(:,:,:)     :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(:,:,:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_igen(VarIn,VarOut,size(VarIn,2)*size(VarIn,3))
+  
+  END SUBROUTINE grid1dTo2d_glo_i2
+  
+  SUBROUTINE grid1dTo2d_glo_i3(VarIn,VarOut)  
+  IMPLICIT NONE  
+    INTEGER,INTENT(IN),DIMENSION(:,:,:,:)     :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(:,:,:,:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_igen(VarIn,VarOut,size(VarIn,2)*size(VarIn,3)*size(VarIn,4))
+  
+  END SUBROUTINE grid1dTo2d_glo_i3
+
+
+  SUBROUTINE grid1dTo2d_glo_r(VarIn,VarOut)  
+  IMPLICIT NONE  
+    REAL,INTENT(IN),DIMENSION(:)     :: VarIn
+    REAL,INTENT(OUT),DIMENSION(:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_rgen(VarIn,VarOut,1)
+  
+  END SUBROUTINE grid1dTo2d_glo_r
+  
+
+  SUBROUTINE grid1dTo2d_glo_r1(VarIn,VarOut)  
+  IMPLICIT NONE  
+    REAL,INTENT(IN),DIMENSION(:,:)     :: VarIn
+    REAL,INTENT(OUT),DIMENSION(:,:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_rgen(VarIn,VarOut,size(VarIn,2))
+  
+  END SUBROUTINE grid1dTo2d_glo_r1
+
+  SUBROUTINE grid1dTo2d_glo_r2(VarIn,VarOut)  
+  IMPLICIT NONE  
+    REAL,INTENT(IN),DIMENSION(:,:,:)     :: VarIn
+    REAL,INTENT(OUT),DIMENSION(:,:,:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_rgen(VarIn,VarOut,size(VarIn,2)*size(VarIn,3))
+  
+  END SUBROUTINE grid1dTo2d_glo_r2
+  
+  SUBROUTINE grid1dTo2d_glo_r3(VarIn,VarOut)  
+  IMPLICIT NONE  
+    REAL,INTENT(IN),DIMENSION(:,:,:,:)     :: VarIn
+    REAL,INTENT(OUT),DIMENSION(:,:,:,:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_rgen(VarIn,VarOut,size(VarIn,2)*size(VarIn,3)*size(VarIn,4))
+  
+  END SUBROUTINE grid1dTo2d_glo_r3
+  
+  
+  
+  SUBROUTINE grid1dTo2d_glo_l(VarIn,VarOut)  
+  IMPLICIT NONE  
+    LOGICAL,INTENT(IN),DIMENSION(:)     :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_lgen(VarIn,VarOut,1)
+  
+  END SUBROUTINE grid1dTo2d_glo_l
+  
+
+  SUBROUTINE grid1dTo2d_glo_l1(VarIn,VarOut)  
+  IMPLICIT NONE  
+    LOGICAL,INTENT(IN),DIMENSION(:,:)     :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(:,:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_lgen(VarIn,VarOut,size(VarIn,2))
+  
+  END SUBROUTINE grid1dTo2d_glo_l1
+
+  SUBROUTINE grid1dTo2d_glo_l2(VarIn,VarOut)  
+  IMPLICIT NONE  
+    LOGICAL,INTENT(IN),DIMENSION(:,:,:)     :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(:,:,:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_lgen(VarIn,VarOut,size(VarIn,2)*size(VarIn,3))
+  
+  END SUBROUTINE grid1dTo2d_glo_l2
+  
+  SUBROUTINE grid1dTo2d_glo_l3(VarIn,VarOut)  
+  IMPLICIT NONE  
+    LOGICAL,INTENT(IN),DIMENSION(:,:,:,:)     :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(:,:,:,:,:)  :: VarOut
+    
+    CALL grid1dTo2d_glo_lgen(VarIn,VarOut,size(VarIn,2)*size(VarIn,3)*size(VarIn,4))
+  
+  END SUBROUTINE grid1dTo2d_glo_l3  
+  
+    SUBROUTINE grid2dTo1d_glo_i(VarIn,VarOut)  
+  IMPLICIT NONE  
+    INTEGER,INTENT(IN),DIMENSION(:,:) :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_igen(VarIn,VarOut,1)
+  
+  END SUBROUTINE grid2dTo1d_glo_i
+  
+
+  SUBROUTINE grid2dTo1d_glo_i1(VarIn,VarOut)  
+  IMPLICIT NONE  
+    INTEGER,INTENT(IN),DIMENSION(:,:,:) :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(:,:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_igen(VarIn,VarOut,size(VarIn,3))
+  
+  END SUBROUTINE grid2dTo1d_glo_i1
+
+  SUBROUTINE grid2dTo1d_glo_i2(VarIn,VarOut)  
+  IMPLICIT NONE  
+    INTEGER,INTENT(IN),DIMENSION(:,:,:,:) :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(:,:,:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_igen(VarIn,VarOut,size(VarIn,3)*size(VarIn,4))
+  
+  END SUBROUTINE grid2dTo1d_glo_i2
+  
+  SUBROUTINE grid2dTo1d_glo_i3(VarIn,VarOut)  
+  IMPLICIT NONE  
+    INTEGER,INTENT(IN),DIMENSION(:,:,:,:,:) :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(:,:,:,:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_igen(VarIn,VarOut,size(VarIn,3)*size(VarIn,4)*size(VarIn,5))
+  
+  END SUBROUTINE grid2dTo1d_glo_i3
+ 
+
+
+
+  SUBROUTINE grid2dTo1d_glo_r(VarIn,VarOut)  
+  IMPLICIT NONE  
+    REAL,INTENT(IN),DIMENSION(:,:) :: VarIn
+    REAL,INTENT(OUT),DIMENSION(:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_rgen(VarIn,VarOut,1)
+  
+  END SUBROUTINE grid2dTo1d_glo_r
+  
+
+  SUBROUTINE grid2dTo1d_glo_r1(VarIn,VarOut)  
+  IMPLICIT NONE  
+    REAL,INTENT(IN),DIMENSION(:,:,:) :: VarIn
+    REAL,INTENT(OUT),DIMENSION(:,:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_rgen(VarIn,VarOut,size(VarIn,3))
+  
+  END SUBROUTINE grid2dTo1d_glo_r1
+
+  SUBROUTINE grid2dTo1d_glo_r2(VarIn,VarOut)  
+  IMPLICIT NONE  
+    REAL,INTENT(IN),DIMENSION(:,:,:,:) :: VarIn
+    REAL,INTENT(OUT),DIMENSION(:,:,:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_rgen(VarIn,VarOut,size(VarIn,3)*size(VarIn,4))
+  
+  END SUBROUTINE grid2dTo1d_glo_r2
+  
+  SUBROUTINE grid2dTo1d_glo_r3(VarIn,VarOut)  
+  IMPLICIT NONE  
+    REAL,INTENT(IN),DIMENSION(:,:,:,:,:) :: VarIn
+    REAL,INTENT(OUT),DIMENSION(:,:,:,:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_rgen(VarIn,VarOut,size(VarIn,3)*size(VarIn,4)*size(VarIn,5))
+  
+  END SUBROUTINE grid2dTo1d_glo_r3
+
+
+
+  SUBROUTINE grid2dTo1d_glo_l(VarIn,VarOut)  
+  IMPLICIT NONE  
+    LOGICAL,INTENT(IN),DIMENSION(:,:) :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_lgen(VarIn,VarOut,1)
+  
+  END SUBROUTINE grid2dTo1d_glo_l
+  
+
+  SUBROUTINE grid2dTo1d_glo_l1(VarIn,VarOut)  
+  IMPLICIT NONE  
+    LOGICAL,INTENT(IN),DIMENSION(:,:,:) :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(:,:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_lgen(VarIn,VarOut,size(VarIn,3))
+  
+  END SUBROUTINE grid2dTo1d_glo_l1
+
+  SUBROUTINE grid2dTo1d_glo_l2(VarIn,VarOut)  
+  IMPLICIT NONE  
+    LOGICAL,INTENT(IN),DIMENSION(:,:,:,:) :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(:,:,:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_lgen(VarIn,VarOut,size(VarIn,3)*size(VarIn,4))
+  
+  END SUBROUTINE grid2dTo1d_glo_l2
+  
+  SUBROUTINE grid2dTo1d_glo_l3(VarIn,VarOut)  
+  IMPLICIT NONE  
+    LOGICAL,INTENT(IN),DIMENSION(:,:,:,:,:) :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(:,:,:,:)  :: VarOut
+    
+    CALL grid2dTo1d_glo_lgen(VarIn,VarOut,size(VarIn,3)*size(VarIn,4)*size(VarIn,5))
+  
+  END SUBROUTINE grid2dTo1d_glo_l3
+
+END MODULE mod_grid_phy_lmdz
+
+
+  
+  SUBROUTINE grid1dTo2d_glo_igen(VarIn,VarOut,dimsize)
+    USE mod_grid_phy_lmdz
+    IMPLICIT NONE
+
+    INTEGER,INTENT(IN) :: dimsize
+    INTEGER,INTENT(IN) ,DIMENSION(klon_glo,dimsize)       :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(nbp_lon*nbp_lat,dimsize)  :: VarOut
+    INTEGER :: i,ij,Offset
+
+    
+    Offset=nbp_lon
+        
+    DO i=1,dimsize
+      DO ij=1,klon_glo
+        VarOut(ij+offset-1,i)=VarIn(ij,i)
+      ENDDO
+    ENDDO
+    
+    
+    DO i=1,dimsize
+      DO ij=1,nbp_lon
+       VarOut(ij,i)=VarIn(1,i)
+      ENDDO
+    ENDDO
+    
+    
+    DO i=1,dimsize
+      DO ij=nbp_lon*(nbp_lat-1)+1,nbp_lat*nbp_lon
+       VarOut(ij,i)=VarIn(klon_glo,i)
+      ENDDO
+    ENDDO
+
+  END SUBROUTINE grid1dTo2d_glo_igen   
+
+
+  SUBROUTINE grid1dTo2d_glo_rgen(VarIn,VarOut,dimsize)
+    USE mod_grid_phy_lmdz
+    IMPLICIT NONE
+
+    INTEGER,INTENT(IN) :: dimsize
+    REAL,INTENT(IN) ,DIMENSION(klon_glo,dimsize)       :: VarIn
+    REAL,INTENT(OUT),DIMENSION(nbp_lon*nbp_lat,dimsize)  :: VarOut
+    INTEGER :: i,ij,Offset
+
+   
+    Offset=nbp_lon
+        
+    DO i=1,dimsize
+      DO ij=1,klon_glo
+        VarOut(ij+offset-1,i)=VarIn(ij,i)
+      ENDDO
+    ENDDO
+    
+    
+    DO i=1,dimsize
+      DO ij=1,nbp_lon
+       VarOut(ij,i)=VarIn(1,i)
+      ENDDO
+    ENDDO
+    
+    
+    DO i=1,dimsize
+      DO ij=nbp_lon*(nbp_lat-1)+1,nbp_lat*nbp_lon
+       VarOut(ij,i)=VarIn(klon_glo,i)
+      ENDDO
+    ENDDO
+
+  END SUBROUTINE grid1dTo2d_glo_rgen   
+
+  SUBROUTINE grid1dTo2d_glo_lgen(VarIn,VarOut,dimsize)
+    USE mod_grid_phy_lmdz
+    IMPLICIT NONE
+    
+    INTEGER,INTENT(IN) :: dimsize
+    LOGICAL,INTENT(IN) ,DIMENSION(klon_glo,dimsize)       :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(nbp_lon*nbp_lat,dimsize)  :: VarOut
+    INTEGER :: i,ij,Offset
+
+    Offset=nbp_lon
+        
+    DO i=1,dimsize
+      DO ij=1,klon_glo
+        VarOut(ij+offset-1,i)=VarIn(ij,i)
+      ENDDO
+    ENDDO
+    
+    
+    DO i=1,dimsize
+      DO ij=1,nbp_lon
+       VarOut(ij,i)=VarIn(1,i)
+      ENDDO
+    ENDDO
+    
+    
+    DO i=1,dimsize
+      DO ij=nbp_lon*(nbp_lat-1)+1,nbp_lat*nbp_lon
+       VarOut(ij,i)=VarIn(klon_glo,i)
+      ENDDO
+    ENDDO
+
+  END SUBROUTINE grid1dTo2d_glo_lgen     
+  
+  
+  SUBROUTINE grid2dTo1d_glo_igen(VarIn,VarOut,dimsize)
+    USE mod_grid_phy_lmdz
+    IMPLICIT NONE
+
+    INTEGER,INTENT(IN) :: dimsize
+    INTEGER,INTENT(IN) ,DIMENSION(nbp_lon*nbp_lat,dimsize) :: VarIn
+    INTEGER,INTENT(OUT),DIMENSION(klon_glo,dimsize)      :: VarOut
+    INTEGER :: i,ij,offset
+
+    offset=nbp_lon
+
+    DO i=1,dimsize
+      DO ij=1,klon_glo
+        VarOut(ij,i)=VarIn(ij+offset-1,i)
+      ENDDO
+    ENDDO
+
+    DO i=1,dimsize
+      VarOut(1,i)=VarIn(1,i)
+    ENDDO
+    
+  END SUBROUTINE grid2dTo1d_glo_igen   
+  
+  SUBROUTINE grid2dTo1d_glo_rgen(VarIn,VarOut,dimsize)
+    USE mod_grid_phy_lmdz
+    IMPLICIT NONE
+
+    INTEGER,INTENT(IN) :: dimsize
+    REAL,INTENT(IN) ,DIMENSION(nbp_lon*nbp_lat,dimsize) :: VarIn
+    REAL,INTENT(OUT),DIMENSION(klon_glo,dimsize)      :: VarOut
+    INTEGER :: i,ij,offset
+
+    offset=nbp_lon
+
+    DO i=1,dimsize
+      DO ij=1,klon_glo
+        VarOut(ij,i)=VarIn(ij+offset-1,i)
+      ENDDO
+    ENDDO
+
+    DO i=1,dimsize
+      VarOut(1,i)=VarIn(1,i)
+    ENDDO
+    
+  END SUBROUTINE grid2dTo1d_glo_rgen 
+    
+  SUBROUTINE grid2dTo1d_glo_lgen(VarIn,VarOut,dimsize)
+    USE mod_grid_phy_lmdz
+    IMPLICIT NONE
+
+    INTEGER,INTENT(IN) :: dimsize
+    LOGICAL,INTENT(IN) ,DIMENSION(nbp_lon*nbp_lat,dimsize) :: VarIn
+    LOGICAL,INTENT(OUT),DIMENSION(klon_glo,dimsize)      :: VarOut
+    INTEGER :: i,ij,offset
+
+    offset=nbp_lon
+
+    DO i=1,dimsize
+      DO ij=1,klon_glo
+        VarOut(ij,i)=VarIn(ij+offset-1,i)
+      ENDDO
+    ENDDO
+
+    DO i=1,dimsize
+      VarOut(1,i)=VarIn(1,i)
+    ENDDO
+    
+  END SUBROUTINE grid2dTo1d_glo_lgen   
Index: /trunk/libf/phytitan/optci.F
===================================================================
--- /trunk/libf/phytitan/optci.F	(revision 101)
+++ /trunk/libf/phytitan/optci.F	(revision 102)
@@ -1,5 +1,5 @@
       SUBROUTINE OPTCI(ykim,nmicro,IPRINT)
+      use dimphy
 #include "dimensions.h"
-#include "dimphy.h"
 #include "microtab.h"
 #include "numchimrad.h"
Index: unk/libf/phytitan/optci.v1
===================================================================
--- /trunk/libf/phytitan/optci.v1	(revision 101)
+++ 	(revision )
@@ -1,549 +1,0 @@
-      SUBROUTINE OPTCI(ykim,nmicro,IPRINT)
-#include "dimensions.h"
-#include "dimphy.h"
-#include "microtab.h"
-#include "numchimrad.h"
-#include "clesphys.h"
-
-c   Arguments:
-c   ---------
-      REAL    ykim(klon,klev,nqmx)
-      integer nmicro
-c   ---------
-
-      PARAMETER(NLAYER=llm,NLEVEL=NLAYER+1)
-      PARAMETER (NSPECI=46,NSPC1I=47,NSPECV=24,NSPC1V=25)
-
-      COMMON /ATM/ Z(NLEVEL),PRESS(NLEVEL),DEN(NLEVEL),TEMP(NLEVEL)
-
-      COMMON /GASS/ CH4(NLEVEL),XN2(NLEVEL),H2(NLEVEL),AR(NLEVEL)
-     & ,XMU(NLEVEL),GAS1(NLAYER),COLDEN(NLAYER)
-
-      COMMON /STRATO/ C2H2(NLAYER),C2H6(NLAYER)
-      COMMON /STRAT2/ HCN(NLAYER)
-
-      COMMON /AERSOL/ RADIUS(NLAYER), XNUMB(NLAYER)
-     & , REALI(NSPECI), XIMGI(NSPECI), REALV(NSPECV), XIMGV(NSPECV)
-
-      COMMON /CLOUD/ RADCLD(NLAYER), XNCLD(NLAYER)
-     & , RCLDI(NSPECI), XICLDI(NSPECI), RCLDV(NSPECV), XICLDV(NSPECV)
-
-      COMMON /TAUS/   TAUHI(klon,NSPECI),TAUCI(klon,NSPECI),
-     &                TAUGI(klon,NSPECI),TAURV(klon,NSPECV),
-     &                TAUHV(klon,NSPECV),TAUCV(klon,NSPECV),
-     &                TAUGV(klon,NSPECV)
-
-      COMMON /TAUD/   TAUHID(klon,NLAYER,NSPECI)
-     &               ,TAUGID(klon,NLAYER,NSPECI)
-     &               ,TAUHVD(klon,NLAYER,NSPECV)
-     &               ,TAUGVD(klon,NLAYER,NSPECV)
-
-
-      COMMON /OPTICI/ DTAUI(klon,NLAYER,NSPECI)
-     &               ,TAUI (klon,NLEVEL,NSPECI)
-     &               ,WBARI(klon,NLAYER,NSPECI)
-     &               ,COSBI(klon,NLAYER,NSPECI)
-
-      COMMON /SPECTI/ BWNI(NSPC1I), WNOI(NSPECI),
-     &                DWNI(NSPECI), WLNI(NSPECI)
-
-      COMMON /PLANT/ CSUBP,RSFI,RSFV,F0PI
-      COMMON /ADJUST/ RHCH4,FH2,FHAZE,FHVIS,FHIR,TAUFAC,RCLOUD,FARGON
-      COMMON /CONST/RGAS,RHOP,PI,SIGMA
-      COMMON /traceurs/qaer
-      COMMON /part/v,r,vrat,dr,dv
-
-      DIMENSION PROD(NLEVEL)
-      real qaer(klon,nlayer,nqmx)
-      real v(nqmx),r(nqmx),vrat,dr(nqmx),dv(nqmx)
-      real xv1(klev,nspeci),xv2(klev,nspeci)
-      real xv3(klev,nspeci)
-      REAL QF1(nqmx,NSPECI),QF2(nqmx,NSPECI)
-      REAL QF3(nqmx,NSPECI),QF4(nqmx,NSPECI)
-      REAL QM1(nqmx,NSPECI),QM2(nqmx,NSPECI)
-      REAL QM3(nqmx,NSPECI),QM4(nqmx,NSPECI)
-      real emu
-      REAL TAEROSM1(NSPECI),TAEROSCATM1(NSPECI),DELTAZM1(NSPECI)
-      
-      save qf1,qf2,qf3,qf4,qm1,qm2,qm3,qm4
-
-      integer iopti,iwarning
-      integer ig_
-      save iopti,iwarning
-      data iopti,iwarning/0,0/
-
-C THE PRESSURE INDUCED TRANSITIONS ARE FROM REGIS
-C THE LAST SEVENTEEN INTERVALS ARE THE BANDS FROM GNF.
-C
-C THIS SUBROUTINE SETS THE OPTICAL CONSTANTS IN THE INFRARED
-C IT CALCUALTES FOR EACH LAYER, FOR EACH SPECRAL INTERVAL IN THE IR
-C LAYER: WBAR, DTAU, COSBAR
-C LEVEL: TAU
-C
-       print*,'START OPTCI'
-c      sum=0.
-c      do nng=2,klon
-c        do i=1,klev
-c         do j=1,nqmx
-c          print*,'j,rj',j,r(j)
-c          print*,'paer',qaer(nng,i,j)
-c           sum=sum+qaer(nng,i,j)*r(j)**3.*1.3333*3.1415*1000.
-c         enddo
-c        enddo
-c        enddo
-c      print*,sum/(klon-1),'SOMME COLONNE/OPTCI'
-
-
-c      do inq=1,nqmx
-c          print*,inq,r(inq),vrat,qaer(12,25,inq)
-c      enddo
-             
-      DO 420 K=1,NSPECI
-C LETS USE THE THOLIN OPTICAL CONSTANTS FOR THE HAZE.
-          CALL THOLIN(WLNI(K),TNR,TNI)
-          REALI(K)=TNR
-          XIMGI(K)=TNI*FHIR
-C SET UP THE OPTICAL CONSTANTS FOR THE CLOUD
-          RCLDI(K)=1.27
-          XICLDI(K)=REFLIQ(WNOI(K))
- 420  CONTINUE
-C
-C
-C ZERO ALL OPTICAL DEPTHS.
-C ??FLAG? FOR SOME APPLCIATIONS THE TOP OPACITY MAY NOT VANISH
-
-c      open  (unit=1,file='xsetupi')
-c      do i=1,klev
-c       read(1,*) a
-c        do j=1,nspeci
-c            read(1,*) xv1(i,j),xv2(i,j),xv3(i,j)
-c        enddo
-c       enddo
-c       close(1)
-
-c************************************************************************
-c************************************************************************
-      DO 79 ig=1,klon      ! BOUCLE SUR GRILLE HORIZONTALE    
-c      print*,'ig NEW optci',ig
-c************************************************************************
-c************************************************************************
-
-        DO 80 K=1,NSPECI
-           TAUHI(ig,K)=0.
-           TAUCI(ig,K)=0.
-           TAUGI(ig,K)=0.
- 80     CONTINUE
-
-        if (.not.ylellouch) then
-        
-            XN2(1) = ykim(ig,1,iradn2)
-            CH4(1) = ykim(ig,1,iradch4)
-             H2(1) = ykim(ig,1,iradh2)
-            do j=2,nlayer
-               XN2(j) = (ykim(ig,j,iradn2)+ykim(ig,j-1,iradn2))/2.
-               CH4(j) = (ykim(ig,j,iradch4)+ykim(ig,j-1,iradch4))/2.
-                H2(j) = (ykim(ig,j,iradh2)+ykim(ig,j-1,iradh2))/2.
-            enddo
-            XN2(nlevel) = ykim(ig,nlayer,iradn2)
-            CH4(nlevel) = ykim(ig,nlayer,iradch4)
-             H2(nlevel) = ykim(ig,nlayer,iradh2)     
-
-            do j=1,nlayer
-               emu = ( xmu(j) + xmu(j+1) )/2.
-               C2H2(j) = ykim(ig,j,iradc2h2) * 26./emu
-               C2H6(j) = ykim(ig,j,iradc2h6) * 30./emu
-                HCN(j) = ykim(ig,j,iradhcn ) * 27./emu
-            enddo
-                  
-        endif
-
-c************************************************************************
-          DO 100 J=1,NLAYER    ! BOUCLE SUR L'ALTITUDE
-c************************************************************************
-
-c      print*,'ig,k,j ',ig,k,j
-
-C SET UP THE COEFFICIENT TO REDUCE MASS PATH TO STP ...SEE NOTES
-C T0 =273.15   PO=1.01325 BAR
-
-        TBAR=0.5*(TEMP(J)+TEMP(J+1))
-        PBAR=SQRT(PRESS(J)*PRESS(J+1))
-        BMU=0.5*(XMU(J+1)+XMU(J))
-        COEF1=RGAS*273.15**2*.5E5* (PRESS(J+1)**2 - PRESS(J)**2)
-     & /(1.01325**2 *EFFG(Z(J))*TBAR*BMU)
-
-      IF (IPRINT .GT. 9) WRITE(6,21) J,EFFG(Z(J)),TBAR,BMU,COEF1
- 21   FORMAT(' J, EFFG, TBAR, BMU, COEF1,: ',I3,1P6E10.3)
-
-c------------------------------------------------------------------------
-         DO 101 K=1,NSPECI     ! BOUCLE SUR LES L.D'O
-c------------------------------------------------------------------------
-
-
-C #1:             HAZE
-C---------------------
-
-C FIRST COMPUTE TAU AEROSOL
-
-
-c
-c                    /\
-c                   /  \
-c                  /    \
-c                 / _O   \
-c                / |/     \
-c               /  / \     \
-c              /   |\ \/\   \
-c             /    || /  \   \
-c             ----------------
-c            |     WARNING    |
-c            |    SLOW DOWN   |
-c             ---------------- 
-
-
-
-
-c*********** EN TRAVAUX ***************************
-
-         TAEROS=0.
-         TAEROSCAT=0.
-         CBAR=0.
-
-      DO inq=1,nmicro         !BOUCLE SUR LES NMICRO TAILLE D"AEROSOLS
-
-
-      IF (WNOI(K).lt.wco) THEN    ! lamda > 56 um
-
-           if (iopti.eq.0) then
-
-c          CALL XMIE(R(inq)*1.e6,REALI(K),XIMGI(K),
-c    &     QEXT,QSCT,QABS,QBAR,WNOI(K))
-
-
-           CALL CMIE(1.E-2/WNOI(K),REALI(K),XIMGI(K),R(inq),
-     &     QEXT,QSCT,QABS,QBAR)
-
-
-           QM1(inq,K)=QEXT
-           QM2(inq,K)=QSCT
-           QM3(inq,K)=QABS
-           QM4(inq,K)=QBAR
-
-          endif         ! end iopti
-
-
-      if (microfi.eq.1) then
-         ig_=ig
-      else
-         ig_=12
-      endif
-
-
-      TAEROS=QM1(inq,K)*qaer(ig_,nlayer+1-J,inq)*1.e-4+TAEROS
-      TAEROSCAT=QM2(inq,K)*qaer(ig_,nlayer+1-J,inq)*1.e-4+TAEROSCAT
-      CBAR=CBAR+QM4(inq,K)*QM2(inq,K)*qaer(ig_,nlayer+1-J,inq)*1.e-4
-
-
-         ELSE                           ! 0.2 < lambda < 56 um
-
-
-            if(R(inq).lt.RF(inq)) THEN
-
-              if (iopti.eq.0) then
-
-                   CALL XMIE(R(inq)*1.e6,REALI(K),XIMGI(K),
-     &             QEXT,QSCT,QABS,QBAR,WNOI(K))
-
-              QM1(inq,K)=QEXT
-              QM2(inq,K)=QSCT
-              QM3(inq,K)=QABS
-              QM4(inq,K)=QBAR
-
-              endif
-
-      if (microfi.eq.1) then
-         ig_=ig
-      else
-         ig_=12
-      endif
-
-
-        TAEROS=QM1(inq,K)*qaer(ig_,nlayer+1-J,inq)*1.e-4+TAEROS
-        TAEROSCAT=QM2(inq,K)*qaer(ig_,nlayer+1-J,inq)*1.e-4+TAEROSCAT
-        CBAR=CBAR+QM4(inq,K)*QM2(inq,K)*qaer(ig_,nlayer+1-J,inq)*1.e-4
-
-           else
-
-               XMONO=(r(inq)/RF(inq))**3.
-               XRULE=1.
-
-            if(XMONO.gt.16384./1.5) then
-             XRULE=(XMONO/16384.)
-             XMONO=16384.
-            endif
-
-
-             if (iopti.eq.0) then
-
-       CALL OPTFRAC(XMONO,10000./WNOI(K)
-     &                         ,QEXT,QSCT,QABS,QBAR)
-
-
-c      CALL CFFFV11(1.e-2/WNOI(K),REALI(K),XIMGI(K),RF(inq),2.
-c    &                ,XMONO,QSCT,QEXT,QABS,QBAR)
-
-
-              QF1(inq,K)=QEXT*XRULE
-              QF2(inq,K)=QSCT*XRULE
-              QF3(inq,K)=QABS*XRULE
-              QF4(inq,K)=QBAR
-
-              endif
-
-      if (microfi.eq.1) then
-         ig_=ig
-      else
-         ig_=12
-      endif
-
-        TAEROS=QF1(inq,K)*qaer(ig_,nlayer+1-J,inq)+TAEROS
-        TAEROSCAT=QF2(inq,K)*qaer(ig_,nlayer+1-J,inq)+TAEROSCAT
-        CBAR=CBAR+QF4(inq,K)*QF2(inq,K)*qaer(ig_,nlayer+1-J,inq)
-
-           endif
-
-               IF(TAEROS.LT.1.e-10) TAEROS=1.e-10  
-
-         ENDIF
-       ENDDO             ! FIN DE LA BOUCLE SUR NMICRO
-
-
-
-
-        CBAR=CBAR/TAEROSCAT
-
-        DELTAZ=Z(J)-Z(J+1)
-
-c --------------------------------------------------------------------
-c profil brume Pascal: fit T (sauf tropopause) et albedo
-c -------------------
-        if( cutoff.eq.1) then
-         IF(PRESS(J).gt.9.e-3) THEN
-          TAEROS=TAEROSM1(K)*DELTAZ/DELTAZM1(K)*0.85
-          TAEROSCAT=TAEROSCATM1(K)*DELTAZ/DELTAZM1(K)*0.85
-c         TAEROS=0.
-c         TAEROSCAT=0.
-         ENDIF
-
-         IF(PRESS(J).gt.1.e-1) THEN
-          TAEROS=TAEROSM1(K)*DELTAZ/DELTAZM1(K)*1.15
-          TAEROSCAT=TAEROSCATM1(K)*DELTAZ/DELTAZM1(K)*1.15
-c         TAEROS=0.
-c         TAEROSCAT=0.
-         ENDIF
-        endif !cutoff=1
-
-c profil brume pour fit T (y compris tropopause), mais ne fit plus albedo...
-c -----------------------
-        if( cutoff.eq.2) then
-         IF(PRESS(J).gt.1.e-1) THEN
-          TAEROS=0.
-          TAEROSCAT=0.
-         ENDIF
-        endif !cutoff=2
-c --------------------------------------------------------------------
-
-         TAEROSM1(K)=TAEROS
-         TAEROSCATM1(K)=TAEROSCAT
-         DELTAZM1(K)=DELTAZ
-
-
-      IF(TAEROSCAT.le.0.) CBAR=0.
-
-c     print*,'HERE, MCKAY AEROSOLS IR'
-c     TAEROS=xv1(j,k)
-c     TAEROSCAT=xv2(j,k)
-c     CBAR=xv3(j,k)
-
-c     if (ig.eq.1) then
-c     if (k.eq.NSPECV/2) then
-c      print*,'@IR',K,J,TAEROS,TAEROSCAT,CBAR
-c     stop'Pour faire des comparaisons'
-c     endif
-c     endif
-
-
-c*********** EN TRAVAUX ***************************
-
-C #2:         CLOUD
-C------------------
-
-C NEXT COMPUTE TAU CLOUD
-      TAUCLD=0.0
-       IF ( XNCLD(J) .GT. 0..and .taufac.gt.0.) THEN
-         print*,'Appel a xmie avec radcld=',radcld(j)
-                CALL XMIE(RADCLD(J),RCLDI(K),XICLDI(K),
-     &                         QEXTC,QSCTC,QABSC,CBARC,WNOI(K))
-                TAUCLD=QEXTC*XNCLD(J)
-        ENDIF
-
-
-C #3:          GAZ
-C------------------
-
-C NOW COMPUTE TAUGAS DUE TO THE PIA TERM ONLY FOR LAMDA LT 940
-       TAUGAS=0.0
-       IF (WNOI(K) .LT. 940. ) THEN
-c        if(ig.eq.1.and.k.eq.nspecv/2) print*,'avant PIA'
-                 CALL PIA(K,TBAR,PNN,PCC,PCN,PHN)
-c        if(ig.eq.1.and.k.eq.nspecv/2) print*,'apres PIA'
-C HERE IS WHERE WE COULD SCALE THE PIA COEFFICEINTS TO FIT DATA
-C BASED ON REGIS' NOTES. ---TGM HAS THIS ADJUST IN IT AS DEFAULT
-                 PCN=PCN*MIN(1.75 , AMAX1(1.0,WNOI(K)/200.))
-C***REPLACE ABOVE WITH: PCN=PCN*1.25*MIN(1.75 , AMAX1(1.0,WNOI(K)/200.))
-C 1.25 FACTOR (NOT FROM DATA) SUGGESTED BY TOON et al. (1988)
-                 TAUGAS=COEF1*
-     &           (XN2(J)*XN2(J)*PNN + CH4(J)*CH4(J)*PCC
-     &           + XN2(J)*CH4(J)*PCN + XN2(J)*H2(J)*PHN)
-            IF (J .EQ. NLAYER .AND. IPRINT .GT. 9)
-     &          WRITE (6,22) WNOI(K),TAUGAS,XN2(J),CH4(J),H2(J),
-     &          TBAR, PNN,PCC,PCN, PHN,
-     &          XN2(J)*XN2(J)*PNN , CH4(J)*CH4(J)*PCC ,
-     &          XN2(J)*CH4(J)*PCN , XN2(J)*H2(J)*PHN
- 22             FORMAT(1X,1P8E10.2)
-       ENDIF
-
-       IF (K .GT. 28) THEN
-                KGAS=K-28
-C     ??FLAG? HERE MUST BE WATCHED CAREFULLY
-                     U=COLDEN(J)*6.02204E23/BMU
-          if(ig.eq.1.and.k.eq.nspecv/2) print*,'Avant GAS2'
-                     if((ylellouch).or.(.not.hcnrad)) then
-                       CALL GAS2_NOHCN(J, KGAS,TBAR,PBAR,U,TAU2)
-                     else
-                       CALL GAS2(J, KGAS,TBAR,PBAR,U,TAU2)
-                     endif
-          if(ig.eq.1.and.k.eq.nspecv/2) print*,'Apres GAS2'
-                     TAUGAS=TAUGAS+TAU2
-       ENDIF
-C
-
-      IF (TAEROS + TAUCLD .GT. 0.) THEN
-         COSBI(ig,J,K)=(CBAR*TAEROSCAT + CBARC*TAUCLD )
-     &                     /(TAEROSCAT + TAUCLD)
-      ELSE
-         COSBI(ig,J,K)=0.0
-      ENDIF
-
-      DTAUI(ig,J,K)=TAUGAS+TAEROS+TAUCLD
-
-      TAUHI(ig,K)=TAUHI(ig,K) + TAEROS
-      TAUHID(ig,J,K)=TAUHI(ig,K)
-
-      TAUGI(ig,K)=TAUGI(ig,K) + TAUGAS
-      TAUGID(ig,J,K)=TAUGI(ig,K)
-
-      TAUCI(ig,K)=TAUCI(ig,K) + TAUCLD
- 
-c     if (ig.eq.1) then
-c     if (k.eq.NSPECI/2.or.k.eq.1.or.k.eq.NSPECI) then
-c      print*,'@IR',K,J,DTAUI(ig,J,K),TAUGAS,TAEROS,TAUCLD         
-c     endif
-c     endif
-
-C ??FLAG? SERIOUS PROBLEM WITH THE CODE HERE!
-
-      TLIMIT=1.E-16
-
-      IF (DTAUI(ig,J,K) .GT. TLIMIT) THEN
-
-c***************** ECHANGE 
-c        WBARI(J,K)=(QSCT*XNUMB(J) + QSCTC*XNCLD(J))
-c****************
-CFC        WBARI(ig,J,K)=(TAEROSCAT + QSCTC*XNCLD(J))
-c****************
-           WBARI(ig,J,K)=TAEROSCAT
-     &   /DTAUI(ig,J,K)
-        if(iwarning.eq.0) 
-     s   print*,'WARNING!!! dans optci xncld non initialise'
-        iwarning=1
-
-      ELSE
-
-         WBARI(ig,J,K)=0.0
-c        PRINT*,'gaz ',taugas,'aerosols',taeros,'nuages',taucld
-c        WRITE (6,999) J,K,DTAUI(ig,J,K)
- 999     FORMAT ('  WARNING TAU MINIMUM AT J,K,DTAU:',2I3,1PE10.3)
-         DTAUI(ig,J,K)=TLIMIT
-
-      ENDIF
-
-c     IF (IPRINT .GT. 9)
-c    & WRITE(6,73)J,K,TAUGAS,TAEROS,QEXT,QSCT
-  73           FORMAT(2I3,1P8E10.3)
- 
-
-c------------------------------------------------------------------------
- 101  CONTINUE   ! FIN BOUCLE L D'O
-c------------------------------------------------------------------------
- 
-      iopti=1
-
-c************************************************************************
- 100  CONTINUE   ! FIN BOUCLE ALTITUDE
-c************************************************************************
- 
-c************************************************************************
-c************************************************************************
-  79  CONTINUE   ! FIN BOUCLE GRILLE HORIZONTALE
-c************************************************************************
-c************************************************************************
-
-
-C TOTAL EXTINCTION OPTICAL DEPTHS
-
-       DO 78 ig=1,klon           ! #1
-
-        DO 119 K=1,NSPECI
-           TAUI(ig,1,K)=0.0
-        DO 119 J=1,NLAYER
-           TAUI(ig,J+1,K)=TAUI(ig,J,K)+DTAUI(ig,J,K)
- 119    CONTINUE
-
-      IF (IPRINT .GT. 2) THEN
-c          WRITE (6,120)
-  120      FORMAT(///'  OPTICAL CONSTANTS IN THE INFRARED')
-
-c        DO 200 K=1,NSPECI           ! #2
-c          WRITE (6,190)
-c          WRITE (6,210)K,WLNI(K),WNOI(K),BWNI(K)
-c    &    ,BWNI(K)+DWNI(K),DWNI(K)
-c          WRITE (6,230)REALI(K),XIMGI(K)
-
-c        DO 195 J=1,NLAYER         !   #3
-c          WRITE (6,220)XNUMB(J), WBARI(ig,J,K),COSBI(ig,J,K)
-c    &                          , DTAUI(ig,J,K),TAUI(ig,J,K)
-c 195    CONTINUE
-
-c        IF(ig.eq.12) WRITE (6,240) TAUI(ig,NLEVEL,K)
-c 200    CONTINUE
-
-          END IF
-  78      CONTINUE
-
-           print*,"TAUI(1400,:,10)=",TAUI(1400,:,10)
-           print*,"DTAUI(1400,:,10)=",DTAUI(1400,:,10)
-
-
-  210 FORMAT(1X,I3,F10.3,F10.2,F10.2,'-',F8.2,F10.3)
-  190 FORMAT(1X//'  SNUM  MICRONS   WAVENU   INTERVAL    DELTA-WN')
-  230 FORMAT(1X,'NREAL(LAYER)= ',1PE10.3,' NIMG(LAYER)= ',E10.3/
-     &' #AEROSOLS   WBAR  COSBAR DTAU TAU')
-  220 FORMAT(5(1X,G9.3))
-  240 FORMAT(41X,G9.3)
-
-c     iopti=1
-      print*, 'FIN OPTCI'
-      stop
-      RETURN
-      END
Index: unk/libf/phytitan/optci.v2
===================================================================
--- /trunk/libf/phytitan/optci.v2	(revision 101)
+++ 	(revision )
@@ -1,233 +1,0 @@
-      SUBROUTINE OPTCI(ykim,nmicro,IPRINT)
-#include "dimensions.h"
-#include "dimphy.h"
-#include "microtab.h"
-#include "numchimrad.h"
-#include "clesphys.h"
-
-c   Arguments:
-c   ---------
-      REAL    ykim(klon,klev,nqmx)
-      integer nmicro
-c   ---------
-
-      PARAMETER(NLAYER=llm,NLEVEL=NLAYER+1)
-      PARAMETER (NSPECI=46,NSPC1I=47,NSPECV=24,NSPC1V=25)
-
-      COMMON /ATM/ Z(NLEVEL),PRESS(NLEVEL),DEN(NLEVEL),TEMP(NLEVEL)
-
-      COMMON /GASS/ CH4(NLEVEL),XN2(NLEVEL),H2(NLEVEL),AR(NLEVEL)
-     & ,XMU(NLEVEL),GAS1(NLAYER),COLDEN(NLAYER)
-
-      COMMON /STRATO/ C2H2(NLAYER),C2H6(NLAYER)
-      COMMON /STRAT2/ HCN(NLAYER)
-
-      COMMON /AERSOL/ RADIUS(NLAYER), XNUMB(NLAYER)
-     & , REALI(NSPECI), XIMGI(NSPECI), REALV(NSPECV), XIMGV(NSPECV)
-
-      COMMON /CLOUD/ RADCLD(NLAYER), XNCLD(NLAYER)
-     & , RCLDI(NSPECI), XICLDI(NSPECI), RCLDV(NSPECV), XICLDV(NSPECV)
-
-      COMMON /TAUS/   TAUHI(klon,NSPECI),TAUCI(klon,NSPECI),
-     &                TAUGI(klon,NSPECI),TAURV(klon,NSPECV),
-     &                TAUHV(klon,NSPECV),TAUCV(klon,NSPECV),
-     &                TAUGV(klon,NSPECV)
-
-      COMMON /TAUD/   TAUHID(klon,NLAYER,NSPECI)
-     &               ,TAUGID(klon,NLAYER,NSPECI)
-     &               ,TAUHVD(klon,NLAYER,NSPECV)
-     &               ,TAUGVD(klon,NLAYER,NSPECV)
-
-
-      COMMON /OPTICI/ DTAUI(klon,NLAYER,NSPECI)
-     &               ,TAUI (klon,NLEVEL,NSPECI)
-     &               ,WBARI(klon,NLAYER,NSPECI)
-     &               ,COSBI(klon,NLAYER,NSPECI)
-
-      COMMON /SPECTI/ BWNI(NSPC1I), WNOI(NSPECI),
-     &                DWNI(NSPECI), WLNI(NSPECI)
-
-      COMMON /PLANT/ CSUBP,RSFI,RSFV,F0PI
-      COMMON /ADJUST/ RHCH4,FH2,FHAZE,FHVIS,FHIR,TAUFAC,RCLOUD,FARGON
-      COMMON /CONST/RGAS,RHOP,PI,SIGMA
-      COMMON /traceurs/qaer
-      COMMON /part/v,r,vrat,dr,dv
-
-      DIMENSION PROD(NLEVEL)
-      real qaer(klon,nlayer,nqmx)
-      real v(nqmx),r(nqmx),vrat,dr(nqmx),dv(nqmx)
-      real xv1(klev,nspeci),xv2(klev,nspeci)
-      real xv3(klev,nspeci)
-      REAL QF1(nqmx,NSPECI),QF2(nqmx,NSPECI)
-      REAL QF3(nqmx,NSPECI),QF4(nqmx,NSPECI)
-      REAL QM1(nqmx,NSPECI),QM2(nqmx,NSPECI)
-      REAL QM3(nqmx,NSPECI),QM4(nqmx,NSPECI)
-      real emu
-      REAL TAEROSM1(NSPECI),TAEROSCATM1(NSPECI),DELTAZM1(NSPECI)
-      
-      save qf1,qf2,qf3,qf4,qm1,qm2,qm3,qm4
-
-      integer iopti,iwarning     ! iopti: premier appel, une seule boucle sur les l.d'o.
-      integer ig_,seulmtunpt
-      save iopti,iwarning,seulmtunpt
-      data iopti,iwarning,seulmtunpt/0,0,0/
-      integer    nmicromax,nmicro1pt
-      parameter (nmicromax=10)
-      real   zqaer_1pt(NLAYER,nmicromax)
-      real   TAUHID_1pt(NLAYER,NSPECI)
-      real   TAUGID_1pt(NLAYER,NSPECI)
-      real   TAUHI_1pt(NSPECI),TAUCI_1pt(NSPECI)
-      real   TAUGI_1pt(NSPECI)
-      real   DTAUI_1pt(NLAYER,NSPECI),TAUI_1pt(NLEVEL,NSPECI)
-      real   WBARI_1pt(NLAYER,NSPECI)
-      real   COSBI_1pt(NLAYER,NSPECI)
-
-C THE PRESSURE INDUCED TRANSITIONS ARE FROM REGIS
-C THE LAST SEVENTEEN INTERVALS ARE THE BANDS FROM GNF.
-C
-C THIS SUBROUTINE SETS THE OPTICAL CONSTANTS IN THE INFRARED
-C IT CALCUALTES FOR EACH LAYER, FOR EACH SPECRAL INTERVAL IN THE IR
-C LAYER: WBAR, DTAU, COSBAR
-C LEVEL: TAU
-C
-       print*,'START OPTCI'
-
-c Diagnostic eventuellement:
-c      if (nmicro.gt.0) then
-c      sum=0.
-c      do nng=2,klon
-c        do i=1,klev
-c         do j=1,nmicro
-c          print*,'j,rj',j,r(j)
-c          print*,'paer',qaer(nng,i,j)
-c           sum=sum+qaer(nng,i,j)*r(j)**3.*1.3333*3.1415*1000.
-c         enddo
-c        enddo
-c        enddo
-c      print*,sum/(klon-1),'SOMME COLONNE/OPTCI'
-c      endif
-
-
-c      do inq=1,nqmx
-c          print*,inq,r(inq),vrat,qaer(12,25,inq)
-c      enddo
-             
-C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-c INITIALISATIONS UNE SEULE FOIS
-C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-
-      if (iopti.eq.0) then
-
-c verif pour taille zqaer_1pt, sachant que si microfi=0 et nqmx=1, 
-c il faut quand meme qu'on lise la look-up table de dim nq=10
-       if (nmicro.gt.nmicromax) then
-          print*,"nmicro.gt.nmicromax",nmicro,nmicromax
-	  print*,"PROBLEME pour zqaer_1pt dans optci !!"
-	  stop
-       endif
-
-      
-      DO 420 K=1,NSPECI
-C LETS USE THE THOLIN OPTICAL CONSTANTS FOR THE HAZE.
-          CALL THOLIN(WLNI(K),TNR,TNI)
-          REALI(K)=TNR
-          XIMGI(K)=TNI*FHIR
-C SET UP THE OPTICAL CONSTANTS FOR THE CLOUD
-          RCLDI(K)=1.27
-          XICLDI(K)=REFLIQ(WNOI(K))
- 420  CONTINUE
-
-C
-C ZERO ALL OPTICAL DEPTHS.
-C ??FLAG? FOR SOME APPLCIATIONS THE TOP OPACITY MAY NOT VANISH
-
-c      open  (unit=1,file='xsetupi')
-c      do i=1,klev
-c       read(1,*) a
-c        do j=1,nspeci
-c            read(1,*) xv1(i,j),xv2(i,j),xv3(i,j)
-c        enddo
-c       enddo
-c       close(1)
-
-      endif    ! fin initialisations premier appel
-
-c************************************************************************
-c************************************************************************
-      DO 79 ig=1,klon      ! BOUCLE SUR GRILLE HORIZONTALE    
-c      print*,'ig NEW optci',ig
-c************************************************************************
-c************************************************************************
-
-        if (.not.ylellouch) then
-        
-            XN2(1) = ykim(ig,1,iradn2)
-            CH4(1) = ykim(ig,1,iradch4)
-             H2(1) = ykim(ig,1,iradh2)
-            do j=2,nlayer
-               XN2(j) = (ykim(ig,j,iradn2)+ykim(ig,j-1,iradn2))/2.
-               CH4(j) = (ykim(ig,j,iradch4)+ykim(ig,j-1,iradch4))/2.
-                H2(j) = (ykim(ig,j,iradh2)+ykim(ig,j-1,iradh2))/2.
-            enddo
-            XN2(nlevel) = ykim(ig,nlayer,iradn2)
-            CH4(nlevel) = ykim(ig,nlayer,iradch4)
-             H2(nlevel) = ykim(ig,nlayer,iradh2)     
-
-            do j=1,nlayer
-               emu = ( xmu(j) + xmu(j+1) )/2.
-               C2H2(j) = ykim(ig,j,iradc2h2) * 26./emu
-               C2H6(j) = ykim(ig,j,iradc2h6) * 30./emu
-                HCN(j) = ykim(ig,j,iradhcn ) * 27./emu
-            enddo
-                  
-        endif
-
-        if (microfi.eq.1) then
-           do iq=1,nmicro
-	      do j=1,NLAYER
-	         zqaer_1pt(j,iq)=qaer(ig,j,iq)
-	      enddo
-	   enddo
-           nmicro1pt=nmicro
-        else
-c initialisation zqaer_1pt a partir d'une look-up table
-c boucle sur nmicromax=10
-            print*,"LOOK-UP TABLE A FAIRE!!  dans optci"
-	    stop
-           nmicro1pt=nmicromax
-        endif
-	
-        if (seulmtunpt.eq.0) then
-	   call optci_1pt(zqaer_1pt,nmicro1pt,iopti,
-     .            COSBI_1pt,DTAUI_1pt,TAUHI_1pt,TAUHID_1pt,TAUCI_1pt,
-     .            TAUGI_1pt,TAUGID_1pt,WBARI_1pt,TAUI_1pt,IPRINT)
-	endif
-
-c Pas de microphysique, ni de composition variable: un seul passage
-c dans optci_1pt.
-        if ((microfi.eq.0).and.(ylellouch)) then
-	   seulmtunpt = 1
-	endif
-	
-	COSBI(ig,:,:)  = COSBI_1pt(:,:) 
-	WBARI(ig,:,:)  = WBARI_1pt(:,:) 
-	DTAUI(ig,:,:)  = DTAUI_1pt(:,:) 
-	TAUHI(ig,:)    = TAUHI_1pt(:) 
-	TAUCI(ig,:)    = TAUCI_1pt(:) 
-	TAUGI(ig,:)    = TAUGI_1pt(:) 
-	TAUI(ig,:,:)   = TAUI_1pt(:,:) 
-	TAUHID(ig,:,:) = TAUHID_1pt(:,:) 
-	TAUGID(ig,:,:) = TAUGID_1pt(:,:) 
-
-c************************************************************************
-c************************************************************************
-  79  CONTINUE   ! FIN BOUCLE GRILLE HORIZONTALE
-c************************************************************************
-c************************************************************************
-
-           print*,"TAUI(1400,:,10)=",TAUI(1400,:,10)
-           print*,"DTAUI(1400,:,10)=",DTAUI(1400,:,10)
-      print*, 'FIN OPTCI'
-      stop
-      RETURN
-      END
Index: /trunk/libf/phytitan/optci_1pt.F
===================================================================
--- /trunk/libf/phytitan/optci_1pt.F	(revision 101)
+++ /trunk/libf/phytitan/optci_1pt.F	(revision 102)
@@ -2,6 +2,6 @@
      .            COSBI_1pt,DTAUI_1pt,TAUHI_1pt,TAUHID_1pt,TAUCI_1pt,
      .            TAUGI_1pt,TAUGID_1pt,WBARI_1pt,TAUI_1pt,IPRINT)
+      use dimphy
 #include "dimensions.h"
-#include "dimphy.h"
 #include "microtab.h"
 #include "numchimrad.h"
Index: /trunk/libf/phytitan/optcv.F
===================================================================
--- /trunk/libf/phytitan/optcv.F	(revision 101)
+++ /trunk/libf/phytitan/optcv.F	(revision 102)
@@ -2,6 +2,6 @@
 
 
+      use dimphy
 #include "dimensions.h"
-#include "dimphy.h"
 #include "microtab.h"
 #include "clesphys.h"
Index: unk/libf/phytitan/optcv.v1
===================================================================
--- /trunk/libf/phytitan/optcv.v1	(revision 101)
+++ 	(revision )
@@ -1,482 +1,0 @@
-      SUBROUTINE OPTCV(nmicro,IPRINT)
-
-
-#include "dimensions.h"
-#include "dimphy.h"
-#include "microtab.h"
-#include "clesphys.h"
-
-c   Argument:
-c   ---------
-      integer nmicro
-c   ---------
-
-      PARAMETER(NLAYER=llm,NLEVEL=NLAYER+1)
-      PARAMETER (NSPECI=46,NSPC1I=47,NSPECV=24,NSPC1V=25)
-      COMMON /ATM/ Z(NLEVEL),PRESS(NLEVEL),DEN(NLEVEL),TEMP(NLEVEL)
-
-      COMMON /GASS/ CH4(NLEVEL),XN2(NLEVEL),H2(NLEVEL),AR(NLEVEL)
-     & ,XMU(NLEVEL),GAS1(NLAYER),COLDEN(NLAYER)
-
-      COMMON /VISGAS/SOLARF(NSPECV),NTERM(NSPECV),PEXPON(NSPECV),
-     &         ATERM(4,NSPECV),BTERM(4,NSPECV)
-
-      COMMON /AERSOL/ RADIUS(NLAYER), XNUMB(NLAYER)
-     & , REALI(NSPECI), XIMGI(NSPECI), REALV(NSPECV), XIMGV(NSPECV)
-
-      COMMON /CLOUD/ RADCLD(NLAYER), XNCLD(NLAYER)
-     &             , RCLDI(NSPECI), XICLDI(NSPECI)
-     &             , RCLDV(NSPECV), XICLDV(NSPECV)
-
-      COMMON /TAUS/   TAUHI(klon,NSPECI), TAUCI(klon,NSPECI)
-     &               ,TAUGI(klon,NSPECI), TAURV(klon,NSPECV)
-     &               ,TAUHV(klon,NSPECV) ,TAUCV(klon,NSPECV)
-     &               ,TAUGV(klon,NSPECV)
-
-      COMMON /TAUD/   TAUHID(klon,NLAYER,NSPECI)
-     &               ,TAUGID(klon,NLAYER,NSPECI)
-     &               ,TAUHVD(klon,NLAYER,NSPECV)
-     &               ,TAUGVD(klon,NLAYER,NSPECV)
-
-      COMMON /OPTICV/ DTAUV(klon,NLAYER,NSPECV,4)
-     &               ,TAUV(klon,NLEVEL,NSPECV,4)
-     &               ,WBARV(klon,NLAYER,NSPECV,4)
-     &               ,COSBV(klon,NLAYER,NSPECV,4)
-
-      COMMON /SPECTV/ BWNV(NSPC1V),WNOV(NSPECV) 
-     &               ,DWNV(NSPECV),WLNV(NSPECV)
-
-      COMMON /PLANT/ CSUBP,RSFI,RSFV,F0PI
-      COMMON /ADJUST/ RHCH4,FH2,FHAZE,FHVIS,FHIR,TAUFAC,RCLOUD,FARGON
-      COMMON /CONST/ RGAS,RHOP,PI,SIGMA
-      COMMON /traceurs/qaer(klon,nlayer,nqmx)
-      COMMON /part/ v(nqmx),r(nqmx),vrat,dr(nqmx),dv(nqmx)
-
-      REAL xv1(klev,NSPECV)
-      REAL xv2(klev,NSPECV)
-      REAL xv3(klev,NSPECV)
-
-      REAL QF1(nqmx,NSPECV),QF2(nqmx,NSPECV)
-      REAL QF3(nqmx,NSPECV),QF4(nqmx,NSPECV)
-      REAL QM1(nqmx,NSPECV),QM2(nqmx,NSPECV)
-      REAL QM3(nqmx,NSPECV),QM4(nqmx,NSPECV)
-
-      save qf1,qf2,qf3,qf4,qm1,qm2,qm3,qm4
- 
-      integer ioptv,iwarning
-      integer ig_
-      save ioptv,iwarning
-      data ioptv,iwarning/0,0/
-
-C*
-C*
-C THIS SUBROUTINE SETS THE OPTICAL CONSTANTS IN THE VISIBLE
-C IT CALCUALTES FOR EACH LAYER, FOR EACH SPECRAL INTERVAL IN THE VIS
-C LAYER: WBAR, DTAU, COSBAR
-C LEVEL: TAU
-C
-       sum=0.
-       PRINT*,'OPTCV'
-       print*,'ATTENTION, TAU UNIFORME DANS OPTCV'
-
-c      do nng=2,klon
-c       do i=1,klev           
-c        do j=1,nqmx
-c          sum=sum+qaer(nng,i,j)*r(j)**3.*1.3333*3.1415*1000.
-c        enddo
-c       enddo
-c       enddo
-c       print*,sum/(klon-1),'SOMME COLONNE/OPTCV'
-
-c      open (unit=1,file='xsetupv')
-c       do j=1,nspecv
-c        read(1,*) a
-c        do i=1,klev
-c            read(1,*) xv1(i,j),xv2(i,j),xv3(i,j)
-c        enddo
-c       enddo
-c       close(1)
-
-      DO 130 K=1,NSPECV
-C LETS USE THE OPTICAL CONSTANTS FOR THOLIN
-      CALL THOLIN(WLNV(K),TNR,TNI)
-      REALV(K)=TNR
-      XIMGV(K)=TNI*FHVIS
-C BUT WE NOW USE THE GEOMETRIC ALBEDO FITTED RESULTS
-C      XIMGV(K)=FITEDT(WLNV(K))
-C      XIMGV(K)=FITEDN(WLNV(K))
-C THE CLOUD IS CLEAR IN THE VISIBLE
-      RCLDV(K)=1.27
-      XICLDV(K)=1.E-7
- 130  CONTINUE
-C
-c******* DEBUT DES BOUCLES ************************
-c     PRINT*, 'AEROSOLS EN VISIBLE'
-
-      DO 101 ig=1,klon       !c! BOUCLE SUR GRILLE HORIZONTALE
-
-      DO 100 K=1,NSPECV         !b! BOUCLE SUR LAMBDA
-
-C ZERO THE COLUMN OPTICAL DEPTHS OF EACH TYPE
-C ??FLAG? THE OPTICAL DEPTH OF THE TOP OF THE MODEL
-C MAY NOT BE ZERO.
-      TAURV(ig,K)=0.
-      TAUHV(ig,K)=0.            ! INTEGRATED TAU.......INITIALIZATION.
-      TAUCV(ig,K)=0.            ! Rayleigh, Haze, Cloud, Gas
-      TAUGV(ig,K)=0.            !   sca,    abs,  abs  , abs
-
-      DO 100 J=1,NLAYER         !a! BOUCLE SUR L"ALTITUDE
-
-C #1:                   HAZE
-c---------------------------
-
-c     CALL THE MIE CODE TO GIVE THE AEROSOL PROPERTIES
-c     USE XFRAC FOR FRACTAL AEROSOLS PROPERTIES AT LAMBDA < 2. um 
-
-
-
-
-c                    /\
-c                   /  \
-c                  /    \
-c                 / _O   \
-c                / |/     \
-c               /  / \     \
-c              /   |\ \/\   \
-c             /    || /  \   \
-c             ----------------
-c            |     WARNING    |
-c            |    SLOW DOWN   |
-c             ---------------- 
-
-
-
-
-c*********** EN TRAVAUX ***************************
- 
-         TAEROS=0.
-         TAEROSCAT=0.
-         CBAR=0.
-
-
-      DO inq=1,nmicro         !BOUCLE SUR LES nmicro TAILLE D"AEROSOLS
-
-
-            IF (R(inq).lt.RF(inq)) THEN    ! aerosols spheriques
-
-              
-            if(ioptv.eq.0.and.J.eq.1) then
-c                  CALL XMIE(R(inq)*1.e6,REALV(K),XIMGV(K),
-c    &             QEXT,QSCT,QABS,QBAR,WNOV(K))
-
-                   CALL CMIE(1.E-2/WNOV(K),REALV(K),XIMGV(K),R(inq),
-     &             QEXT,QSCT,QABS,QBAR)
-
-c       print*,'inq=',inq,' QM1=',QM1(inq,K),' QEXT=',QEXT
-
-              QM1(inq,K)=QEXT
-              QM2(inq,K)=QSCT
-              QM3(inq,K)=QABS
-              QM4(inq,K)=QBAR
-            endif
-
-
-      if (microfi.eq.1) then
-         ig_=ig
-      else
-         ig_=12
-      endif
-
-       TAEROS=QM1(inq,K)*qaer(ig_,NLAYER+1-J,inq)*1.e-4+TAEROS
-       TAEROSCAT=QM2(inq,K)*qaer(ig_,NLAYER+1-J,inq)*1.e-4+TAEROSCAT
-       CBAR=CBAR+QM4(inq,K)*QM2(inq,K)*qaer(ig_,NLAYER+1-J,inq)*1.e-4
-
-            ELSE                        ! aerosols fractals
-
-               XMONO=(R(inq)/RF(inq))**3.
-               XRULE=1. 
-
-            if(XMONO.gt.16384./1.5) then 
-             XRULE=(XMONO/16384.) 
-             XMONO=16384.
-            endif 
-
-            if(ioptv.eq.0.and.J.eq.1) then
-
-        CALL OPTFRAC(XMONO,10000./WNOV(K)
-     &                        ,QEXT,QSCT,QABS,QBAR)
-
-c       CALL CFFFV11(1.e-2/WNOV(K),REALV(K),XIMGV(K),RF(inq),2.
-c    &   ,XMONO,QSCT,QEXT,QABS,QBAR)
-
-
-              QF1(inq,K)=QEXT*XRULE
-              QF2(inq,K)=QSCT*XRULE
-              QF3(inq,K)=QABS*XRULE
-              QF4(inq,K)=QBAR
-
-c       print*,'inq=',inq,' QF1=',QF1(inq,K),' QEXT=',QEXT,' XRULE=',XRULE
-                   
-            endif
-
-
-      if (microfi.eq.1) then
-         ig_=ig
-      else
-         ig_=12
-      endif
-
-        TAEROS=QF1(inq,K)*qaer(ig_,NLAYER+1-J,inq)+TAEROS
-        TAEROSCAT=QF2(inq,K)*qaer(ig_,NLAYER+1-J,inq)+TAEROSCAT
-        CBAR=CBAR+QF4(inq,K)*QF2(inq,K)*qaer(ig_,NLAYER+1-J,inq)
-
-
-           ENDIF 
-
-
-       ENDDO    ! nmicro
-
-
-       CBAR=CBAR/TAEROSCAT
-
-        DELTAZ=Z(J)-Z(J+1)
-
-c --------------------------------------------------------------------
-c profil brume Pascal: fit T (sauf tropopause) et albedo
-c -------------------
-        if( cutoff.eq.1) then
-         IF(PRESS(J).gt.9.e-3) THEN
-          TAEROS=TAEROSM1*DELTAZ/DELTAZM1*0.85
-          TAEROSCAT=TAEROSCATM1*DELTAZ/DELTAZM1*0.85
-c         TAEROS=0.
-c         TAEROSCAT=0.
-         ENDIF
-
-         IF(PRESS(J).gt.1.e-1) THEN
-          TAEROS=TAEROSM1*DELTAZ/DELTAZM1*1.15
-          TAEROSCAT=TAEROSCATM1*DELTAZ/DELTAZM1*1.15
-c         TAEROS=0.
-c         TAEROSCAT=0.
-         ENDIF
-        endif !cutoff=1
-
-c profil brume pour fit T (y compris tropopause), mais ne fit plus albedo...
-c -----------------------
-        if( cutoff.eq.2) then
-         IF(PRESS(J).gt.1.e-1) THEN
-          TAEROS=0.
-          TAEROSCAT=0.
-         ENDIF
-        endif !cutoff=2
-c --------------------------------------------------------------------
-
-         TAEROSM1=TAEROS
-         TAEROSCATM1=TAEROSCAT
-         DELTAZM1=DELTAZ
-
-
-       IF (TAEROSCAT.le.0.) CBAR=0. 
-
-c          if(ig.eq.12) then
-c          if(j.eq.1) print*,'NEWK',wlnv(k)
-c          print*,j,TAEROS,xv1(j,k),'   ', TAEROSCAT/TAEROS,
-c    &     xv2(j,k)/xv1(j,k),'  ',CBAR,xv3(j,k)
-c          print*,' ' 
-c          endif
-           
-
-c      print*,'HERE; MCKAY AEROSOLS VIS'
-
-c      TAEROSCAT=xv2(j,k)
-c      TAEROS=xv1(j,k)
-c      CBAR=xv3(j,k)
-
-c     if (ig.eq.1) then
-c     if (k.eq.NSPECV/2) then   
-c      print*,'@VI',K,J,TAEROS,TAEROSCAT,CBAR
-c      print*,'@  ',K,J,QF1(1,K),QF2(1,K),qaer(12,NLAYER+1-J,1)
-c      print*,'@  ',K,J,QF1(3,K),QF2(3,K),qaer(12,NLAYER+1-J,3)
-c      print*,'@  ',K,J,QF1(5,K),QF2(5,K),qaer(12,NLAYER+1-J,5)
-c      print*,'@  ',K,J,QF1(7,K),QF2(7,K),qaer(12,NLAYER+1-J,7)
-c      print*,'@  ',K,J,QF1(9,K),QF2(9,K),qaer(12,NLAYER+1-J,9)
-c      print*
-c     endif
-c     endif
-
-
-
-c*********** EN TRAVAUX ***************************
-
-C #2:                   RAYLEIGH
-c-------------------------------
-
-C RAYLEIGH SCATTERING STRAIGHT FROM HANSEN AND TRAVIS...SEE NOTES
-C RATIOED BY THE LAYER COLUMN NUMBER TO THE TOTAL
-C COLUMN NUMBER ON EARTH. CM-2
-C THIS IS THE SCATTERING BY THE ATMOSPHERE
-
-      TAURAY=(COLDEN(J)*28.9/(XMU(J)*1013.25))*
-     &(.008569/WLNV(K)**4)*(1.+.0113/WLNV(K)**2+.00013/WLNV(K)**4)
-
-c       PRINT*,WLNV(K)
-c      COLX=0.
-c      COLP=0.
-c      COLT=0.
-c     DO IU=1,NLAYER
-c      COLP=COLDEN(IU)*1.e+1*1.35+COLP
-c     TAURAY=(COLDEN(IU)*28.9/(XMU(IU)*1013.25))*
-c    & (.008569/WLNV(K)**4)*(1.+.0113/WLNV(K)**2
-c    & +.00013/WLNV(K)**4)
-c      COLT=COLT+TAURAY
-c      COLX=COLDEN(IU)*1.e+1/(1.E5*28./22.4E3)*1.e-1*0.0933e-1+COLX
-c                           |    
-c                           |    
-c           g/cm2->kg/m2    |  m2/kg   
-c      Print*,IU, tauray, 
-c    &   COLDEN(IU)*1.e+1/(1.E5*28./22.4E3)*1.e-1*0.543e-1
-c     ENDDO
-c       PRINT*,COLP,' PRESSURE AT GROUND;'
-c       PRINT*,COLX,' TAU_GAS AT GROUND;'
-c       print*,colt,colx,' COLT, COLX'
-c      STOP
-							
-c       DZ=Z(J)-Z(J+1)
-c     PRINT*, Z(J),WLNV(K),
-c    &(28.9/(XMU(J)*1013.25))*(.008569/WLNV(K)**4)*
-c    &(1.+.0113/WLNV(K)**2+.00013/WLNV(K)**4)
-c    & ,COLDEN(J)/DZ/100000.,
-c    &(28.9/(XMU(J)*1013.25))*(.008569/WLNV(K)**4)*
-c    &(1.+.0113/WLNV(K)**2+.00013/WLNV(K)**4)
-c    & *COLDEN(J)/DZ/100000.
-    
- 
-
-C #3:                   CLOUD
-c----------------------------
-
-C NEXT COMPUTE TAU CLOUD
-
-      TAUCLD=0.0
-c             XNCLD(J)=0.
-      IF ( XNCLD(J) .GT. 0. .and .taufac.gt.0.) THEN
-                CALL XMIE(RADCLD(J),RCLDV(K),XICLDV(K),
-     &                         QEXTC,QSCTC,QABSC,CBARC,WNOV(K))
-                TAUCLD=QEXTC*XNCLD(J)         
-      ENDIF
-C
-      TAURV(ig,K)=TAURV(ig,K)+TAURAY
-      TAUGVD(ig,J,K)=TAURV(ig,K)
-
-      TAUHV(ig,K)=TAUHV(ig,K)+TAEROS          ! INTEGRATED Quant.
-      TAUHVD(ig,J,K)=TAUHV(ig,K)
-
-      TAUCV(ig,K)=TAUCV(ig,K)+TAUCLD
-
-C #4:                  TAUGAS
-C----------------------------
-
-C LOOP OVER THE NTERMS
-C THIS IS THE ABSORPTION BY THE ATMOSPHERE (METHANE)
-
-
-      DO 909 NT=1,NTERM(K)
-      TAUGAS=COLDEN(J)*GAS1(J)*BTERM(NT,K)*
-     &  (   (PRESS(J+1) + PRESS(J))*.5  )**PEXPON(K)
-
-
-C COMPUTE THE AVERAGE COSBAR AND WBAR
-C&&
-
-c     CBAR=MIN(1.0,1.05*CBAR)       ! THE HAZE FORWARD SCATTERING 5%(WHY?)
-      COSBV(ig,J,K,NT)=(CBAR*TAEROSCAT + CBARC*TAUCLD)
-     &  /(TAEROSCAT+TAUCLD+TAURAY)     !CBAR_RAY=0.
-c        print*,'CBV',J,K,NT,CBAR,TAEROSCAT,CBARC,TAUCLD       
-
-      DTAUV(ig,J,K,NT)=TAUGAS+TAEROS+TAURAY+TAUCLD       !TOTAL TAU_EXT
-      TAUGV(ig,K)=TAUGV(ig,K)+TAUGAS*ATERM(NT,K)         !TAU_ABS_METH INTEG.
-
-C WE LET W RAYLEIGH BE .999 OR W=1 WHEN ONLY RAYLEIGH=PROBLEM FOR TRID
-c WE HAVE ASSUMED ABOVE THAT COSBAR FOR RAYLEIGH IS ZERO.
-c     if (ig.eq.1) then
-c     if (k.eq.NSPECV/2) then   
-c      print*,'@VI',K,J,DTAUV(ig,J,K,1),TAUGAS,TAEROS,TAUCLD
-c     endif
-c     endif
-
-
-c***************** ECHANGE 
-c     WBARV(J,K,NT)=(QSCT*XNUMB(J)+TAURAY*0.9999999 + QSCTC*XNCLD(J) )
-c****************
-      WBARV(ig,J,K,NT)=(TAEROSCAT+TAURAY*0.9999999 + QSCTC*XNCLD(J) )
-c     WBARV(ig,J,K,NT)=(TAEROSCAT+TAURAY*0.9999999 )
-     &            /(TAUGAS+TAEROS+TAURAY+TAUCLD)
-c**************** 
-      IF((TAEROS+TAUCLD+TAURAY+TAUCLD).le.0.) WBARV(ig,J,K,NT)=0.
-      IF((TAEROS+TAUCLD+TAURAY).le.0.) COSBV(ig,J,K,NT)=0.
-
-c     print*,'WBV',J,K,NT,TAEROSCAT,TAURAY,QSCTC*XNCLD(J)
-c     print*,'WBV',J,K,NT,TAEROS,TAUGAS,TAURAY,TAUCLD
-c     print*,Z(j),J,K,NT,TAUV(1,j,K,NT),WBARV(1,j,K,NT),COSBV(1,j,K,NT)
-
- 909  CONTINUE
-      TAUGVD(ig,J,K)=TAUGVD(ig,J,K)+TAUGV(ig,K)
- 100  CONTINUE
-       ioptv=1
- 101  CONTINUE
-c HERE END OF THE LOOPS *******
-c******************************
-         
-        DO 102 ig=1,klon
-
-C TOTAL EXTINCTION OPTICAL DEPTHS
-          DO 119 K=1,NSPECV
-C LOOP OVER NTERMS
-           DO 119 NT=1,NTERM(K)
-           TAUV(ig,1,K,NT)=0.0
-             DO 119 J=1,NLAYER
-             TAUV(ig,J+1,K,NT)=TAUV(ig,J,K,NT)+DTAUV(ig,J,K,NT)
- 119     CONTINUE
-
-c       print*,'SETUP'
-c      do i=1,NSPECV
-c      print*,WLNV(i)
-c       do j=1,NLAYER+1
-c       print*,Z(j),TAUV(1,j,i,1),WBARV(1,j,i,1),COSBV(1,j,i,1)
-c       enddo
-c      enddo
-c
-c     IF (IPRINT .GT. 1) THEN
-c           NT=1
-c     IF (2 .GT. 1) THEN
-c          WRITE (6,120)
-c 120      FORMAT(///'  OPTICAL CONSTANTS IN THE VISIBLE (@EQUATOR) ')
-c          WRITE(6,*) 'latitude:',ig
-c          DO 200 K=1,NSPECV
-c          WRITE (6,190)
-c          WRITE (6,210)K,WLNV(K),WNOV(K),BWNV(K)
-c    &    ,BWNV(K)+DWNV(K),DWNV(K)
-c          WRITE (6,230)REALV(K),XIMGV(K)
-c          DO 195 J=1,NLAYER,NLAYER
-C RECALCULATE FOR PRINT OUT ONLY, ONLY FIRST NTERM AT ig=12 (EQUATOR)
-c          WRITE (6,220)XNUMB(J), WBARV(ig,J,K,NT),COSBV(ig,J,K,NT)
-c    &      ,DTAUV(ig,J,K,NT),TAUV(ig,J,K,NT)
-c 195      CONTINUE
-c          WRITE (6,240) TAUV(IG,NLEVEL,K,NT)
-c 200      CONTINUE
-c     END IF
-  102 CONTINUE
-  210 FORMAT(1X,I3,F10.3,F10.2,F10.2,'-',F8.2,F10.3)
-  190 FORMAT(1X//'  SNUM  MICRONS   WAVENU   INTERVAL    DELTA-WN')
-  230 FORMAT(1X,'NREAL(LAYER)= ',1PE10.3,' NIMG(LAYER)= ',E10.3/
-     &' #AEROSOLS   WBAR  COSBAR       DTAU     TAU'
-     & ,9X,'RAY     GAS    AEROSOL')
-  220 FORMAT(8(1X,F9.3))
-  240 FORMAT(41X,F9.3)
-
-           print*,"TAUV(1400,:,10,2)=",TAUV(1400,:,10,2)
-           print*,"DTAUV(1400,:,10,2)=",DTAUV(1400,:,10,2)
-c      ioptv=1
-       PRINT*, 'FIN OPTCV'
-       stop
-      RETURN
-      END
Index: unk/libf/phytitan/optcv.v2
===================================================================
--- /trunk/libf/phytitan/optcv.v2	(revision 101)
+++ 	(revision )
@@ -1,193 +1,0 @@
-      SUBROUTINE OPTCV(nmicro,IPRINT)
-
-
-#include "dimensions.h"
-#include "dimphy.h"
-#include "microtab.h"
-#include "clesphys.h"
-
-c   Argument:
-c   ---------
-      integer nmicro
-c   ---------
-
-      PARAMETER(NLAYER=llm,NLEVEL=NLAYER+1)
-      PARAMETER (NSPECI=46,NSPC1I=47,NSPECV=24,NSPC1V=25)
-      COMMON /ATM/ Z(NLEVEL),PRESS(NLEVEL),DEN(NLEVEL),TEMP(NLEVEL)
-
-      COMMON /GASS/ CH4(NLEVEL),XN2(NLEVEL),H2(NLEVEL),AR(NLEVEL)
-     & ,XMU(NLEVEL),GAS1(NLAYER),COLDEN(NLAYER)
-
-      COMMON /VISGAS/SOLARF(NSPECV),NTERM(NSPECV),PEXPON(NSPECV),
-     &         ATERM(4,NSPECV),BTERM(4,NSPECV)
-
-      COMMON /AERSOL/ RADIUS(NLAYER), XNUMB(NLAYER)
-     & , REALI(NSPECI), XIMGI(NSPECI), REALV(NSPECV), XIMGV(NSPECV)
-
-      COMMON /CLOUD/ RADCLD(NLAYER), XNCLD(NLAYER)
-     &             , RCLDI(NSPECI), XICLDI(NSPECI)
-     &             , RCLDV(NSPECV), XICLDV(NSPECV)
-
-      COMMON /TAUS/   TAUHI(klon,NSPECI), TAUCI(klon,NSPECI)
-     &               ,TAUGI(klon,NSPECI), TAURV(klon,NSPECV)
-     &               ,TAUHV(klon,NSPECV) ,TAUCV(klon,NSPECV)
-     &               ,TAUGV(klon,NSPECV)
-
-      COMMON /TAUD/   TAUHID(klon,NLAYER,NSPECI)
-     &               ,TAUGID(klon,NLAYER,NSPECI)
-     &               ,TAUHVD(klon,NLAYER,NSPECV)
-     &               ,TAUGVD(klon,NLAYER,NSPECV)
-
-      COMMON /OPTICV/ DTAUV(klon,NLAYER,NSPECV,4)
-     &               ,TAUV(klon,NLEVEL,NSPECV,4)
-     &               ,WBARV(klon,NLAYER,NSPECV,4)
-     &               ,COSBV(klon,NLAYER,NSPECV,4)
-
-      COMMON /SPECTV/ BWNV(NSPC1V),WNOV(NSPECV) 
-     &               ,DWNV(NSPECV),WLNV(NSPECV)
-
-      COMMON /PLANT/ CSUBP,RSFI,RSFV,F0PI
-      COMMON /ADJUST/ RHCH4,FH2,FHAZE,FHVIS,FHIR,TAUFAC,RCLOUD,FARGON
-      COMMON /CONST/ RGAS,RHOP,PI,SIGMA
-      COMMON /traceurs/qaer(klon,nlayer,nqmx)
-      COMMON /part/ v(nqmx),r(nqmx),vrat,dr(nqmx),dv(nqmx)
-
-      REAL xv1(klev,NSPECV)
-      REAL xv2(klev,NSPECV)
-      REAL xv3(klev,NSPECV)
-
-      REAL QF1(nqmx,NSPECV),QF2(nqmx,NSPECV)
-      REAL QF3(nqmx,NSPECV),QF4(nqmx,NSPECV)
-      REAL QM1(nqmx,NSPECV),QM2(nqmx,NSPECV)
-      REAL QM3(nqmx,NSPECV),QM4(nqmx,NSPECV)
-
-      save qf1,qf2,qf3,qf4,qm1,qm2,qm3,qm4
- 
-      integer ioptv,iwarning     ! ioptv: premier appel, une seule boucle sur les l.d'o.
-      integer ig_,seulmtunpt
-      save ioptv,iwarning,seulmtunpt
-      data ioptv,iwarning,seulmtunpt/0,0,0/
-      integer    nmicromax,nmicro1pt
-      parameter (nmicromax=10)
-      real   zqaer_1pt(NLAYER,nmicromax)
-      real   TAUHVD_1pt(NLAYER,NSPECV)
-      real   TAUGVD_1pt(NLAYER,NSPECV)
-      real   TAUHV_1pt(NSPECV),TAUCV_1pt(NSPECV)
-      real   TAURV_1pt(NSPECV),TAUGV_1pt(NSPECV)
-      real   DTAUV_1pt(NLAYER,NSPECV,4),TAUV_1pt(NLEVEL,NSPECV,4)
-      real   WBARV_1pt(NLAYER,NSPECV,4)
-      real   COSBV_1pt(NLAYER,NSPECV,4)
-
-C*
-C THIS SUBROUTINE SETS THE OPTICAL CONSTANTS IN THE VISIBLE
-C IT CALCUALTES FOR EACH LAYER, FOR EACH SPECRAL INTERVAL IN THE VIS
-C LAYER: WBAR, DTAU, COSBAR
-C LEVEL: TAU
-C
-       sum=0.
-       PRINT*,'OPTCV'
-       print*,'ATTENTION, TAU UNIFORME DANS OPTCV'
-
-c      do nng=2,klon
-c       do i=1,klev           
-c        do j=1,nqmx
-c          sum=sum+qaer(nng,i,j)*r(j)**3.*1.3333*3.1415*1000.
-c        enddo
-c       enddo
-c       enddo
-c       print*,sum/(klon-1),'SOMME COLONNE/OPTCV'
-
-             
-C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-c INITIALISATIONS UNE SEULE FOIS
-C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-
-      if (ioptv.eq.0) then
-
-c verif pour taille zqaer_1pt, sachant que si microfi=0 et nqmx=1, 
-c il faut quand meme qu'on lise la look-up table de dim nq=10
-       if (nmicro.gt.nmicromax) then
-          print*,"nmicro.gt.nmicromax",nmicro,nmicromax
-	  print*,"PROBLEME pour zqaer_1pt dans optcv !!"
-	  stop
-       endif
-
-      
-      DO 130 K=1,NSPECV
-C LETS USE THE OPTICAL CONSTANTS FOR THOLIN
-      CALL THOLIN(WLNV(K),TNR,TNI)
-      REALV(K)=TNR
-      XIMGV(K)=TNI*FHVIS
-C BUT WE NOW USE THE GEOMETRIC ALBEDO FITTED RESULTS
-C      XIMGV(K)=FITEDT(WLNV(K))
-C      XIMGV(K)=FITEDN(WLNV(K))
-C THE CLOUD IS CLEAR IN THE VISIBLE
-      RCLDV(K)=1.27
-      XICLDV(K)=1.E-7
- 130  CONTINUE
-C
-c      open (unit=1,file='xsetupv')
-c       do j=1,nspecv
-c        read(1,*) a
-c        do i=1,klev
-c            read(1,*) xv1(i,j),xv2(i,j),xv3(i,j)
-c        enddo
-c       enddo
-c       close(1)
-
-      endif    ! fin initialisations premier appel
-
-c******* DEBUT DES BOUCLE GRILLE ************************
-c     PRINT*, 'AEROSOLS EN VISIBLE'
-
-      DO 101 ig=1,klon       !c! BOUCLE SUR GRILLE HORIZONTALE
-
-        if (microfi.eq.1) then
-           do iq=1,nmicro
-	      do j=1,NLAYER
-	         zqaer_1pt(j,iq)=qaer(ig,j,iq)
-	      enddo
-	   enddo
-           nmicro1pt=nmicro
-        else
-c initialisation zqaer_1pt a partir d'une look-up table
-c boucle sur nmicromax=10
-            print*,"LOOK-UP TABLE A FAIRE!!  dans optcv"
-	    stop
-           nmicro1pt=nmicromax
-        endif
-	
-        if (seulmtunpt.eq.0) then
-	   call optcv_1pt(zqaer_1pt,nmicro1pt,ioptv,
-     .            COSBV_1pt,DTAUV_1pt,TAUHV_1pt,TAUHVD_1pt,TAUCV_1pt,
-     .       TAURV_1pt,TAUGV_1pt,TAUGVD_1pt,WBARV_1pt,TAUV_1pt,IPRINT)
-	endif
-
-c Pas de microphysique, ni de composition variable: un seul passage
-c dans optci_1pt.
-        if ((microfi.eq.0).and.(ylellouch)) then
-	   seulmtunpt = 1
-	endif
-	
-	COSBV(ig,:,:,:)= COSBV_1pt(:,:,:) 
-	WBARV(ig,:,:,:)= WBARV_1pt(:,:,:) 
-	DTAUV(ig,:,:,:)= DTAUV_1pt(:,:,:) 
-	TAUHV(ig,:)    = TAUHV_1pt(:) 
-	TAUCV(ig,:)    = TAUCV_1pt(:) 
-	TAURV(ig,:)    = TAURV_1pt(:) 
-	TAUGV(ig,:)    = TAUGV_1pt(:) 
-	TAUV(ig,:,:,:) = TAUV_1pt(:,:,:) 
-	TAUHVD(ig,:,:) = TAUHVD_1pt(:,:) 
-	TAUGVD(ig,:,:) = TAUGVD_1pt(:,:) 
-
- 101  CONTINUE
-
-c FIN BOUCLE GRILLE     *******
-c******************************
-         
-           print*,"TAUV(1400,:,10,2)=",TAUV(1400,:,10,2)
-           print*,"DTAUV(1400,:,10,2)=",DTAUV(1400,:,10,2)
-       PRINT*, 'FIN OPTCV'
-       stop
-      RETURN
-      END
Index: /trunk/libf/phytitan/optcv_1pt.F
===================================================================
--- /trunk/libf/phytitan/optcv_1pt.F	(revision 101)
+++ /trunk/libf/phytitan/optcv_1pt.F	(revision 102)
@@ -4,6 +4,6 @@
 
 
+      use dimphy
 #include "dimensions.h"
-#include "dimphy.h"
 #include "microtab.h"
 #include "clesphys.h"
Index: /trunk/libf/phytitan/pg2.F
===================================================================
--- /trunk/libf/phytitan/pg2.F	(revision 101)
+++ /trunk/libf/phytitan/pg2.F	(revision 102)
@@ -41,11 +41,11 @@
 c------------------------------------------------------
 
+      use dimphy
+      USE comgeomphy
 #include "dimensions.h"
-#include "dimphy.h"
 #include "microtab.h"
 #include "paramet.h"
 #include "aerprod.h"
 #include "clesphys.h"
-#include "comgeomphy.h"
 #include "YOMCST.h"
 
Index: /trunk/libf/phytitan/phyetat0.F
===================================================================
--- /trunk/libf/phytitan/phyetat0.F	(revision 101)
+++ /trunk/libf/phytitan/phyetat0.F	(revision 102)
@@ -10,11 +10,11 @@
      .           tabcntr0,
      .           t_ancien,ancien_ok)
-      IMPLICIT none
 c======================================================================
 c Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
 c Objet: Lecture de l'etat initial pour la physique
 c======================================================================
+      use dimphy
+      IMPLICIT none
 #include "dimensions.h"
-#include "dimphy.h"
 #include "netcdf.inc"
 #include "dimsoil.h"
Index: /trunk/libf/phytitan/phyredem.F
===================================================================
--- /trunk/libf/phytitan/phyredem.F	(revision 101)
+++ /trunk/libf/phytitan/phyredem.F	(revision 102)
@@ -11,4 +11,5 @@
 
       USE control_mod
+      use dimphy
       IMPLICIT none
 c======================================================================
@@ -17,5 +18,4 @@
 c======================================================================
 #include "dimensions.h"
-#include "dimphy.h"
 #include "netcdf.inc"
 #include "dimsoil.h"
@@ -75,7 +75,7 @@
       IF( cycle_diurne ) tab_cntrl( 7 ) = 1.
       IF(   soil_model ) tab_cntrl( 8 ) = 1.
-      IF(     new_oliq ) tab_cntrl( 9 ) = 1.
       IF(     ok_orodr ) tab_cntrl(10 ) = 1.
       IF(     ok_orolf ) tab_cntrl(11 ) = 1.
+      IF( ok_gw_nonoro ) tab_cntrl(12 ) = 1.
 
       tab_cntrl(13) = day_end
Index: /trunk/libf/phytitan/physiq.F
===================================================================
--- /trunk/libf/phytitan/physiq.F	(revision 101)
+++ /trunk/libf/phytitan/physiq.F	(revision 102)
@@ -10,9 +10,4 @@
      .            d_u, d_v, d_t, d_qx, d_ps)
 
-      USE ioipsl
-      USE histcom
-      USE infotrac
-      USE control_mod
-      IMPLICIT none
 c======================================================================
 c
@@ -29,11 +24,4 @@
 cAA                  -  stockage des moyennes des champs necessaires
 cAA                     en mode traceur off-line 
-c======================================================================
-c   CLEFS CPP POUR LES IO
-c   =====================
-#define histmth
-#define histday
-#define histins
-c======================================================================
 c    modif   ( P. Le Van ,  12/10/98 )
 c
@@ -67,14 +55,25 @@
 c d_ps----output-R-tendance physique de la pression au sol
 c======================================================================
+      USE ioipsl
+      USE histcom
+      USE infotrac
+      USE control_mod
+      use dimphy
+      USE comgeomphy
+      IMPLICIT none
+c======================================================================
+c   CLEFS CPP POUR LES IO
+c   =====================
+#define histmth
+#define histday
+#define histins
+c======================================================================
 #include "dimensions.h"
       integer jjmp1
       parameter (jjmp1=jjm+1-1/jjm)
-#include "dimphy.h"
 #include "dimsoil.h"
 #include "clesphys.h"
 #include "temps.h"
-#include "comgeomphy.h"
 #include "iniprint.h"
-#include "raddim.h"
 #include "timerad.h" 
 #include "logic.h"
@@ -583,5 +582,5 @@
       do iq=1,nqmax
          nom = tname(iq)
-c        print*,iq,"nom=",nom,"tname=",tname(iq),"tnom=",tnom(iq)
+c        print*,iq,"nom=",nom,"tname=",tname(iq)
          print*,iq,"nom=",nom
          if (nom(1:1).eq."q") then
Index: /trunk/libf/phytitan/phytrac.F
===================================================================
--- /trunk/libf/phytitan/phytrac.F	(revision 101)
+++ /trunk/libf/phytitan/phytrac.F	(revision 102)
@@ -3,8 +3,4 @@
      .                   pplev,pplay,delp,ptemp,pmu0,pfract,pdecli,
      .                   lonsol,tr_seri,d_tr_mph,d_tr_kim)
-
-      USE infotrac
-
-      IMPLICIT none
 
 c======================================================================
@@ -32,6 +28,8 @@
 c d_tr_kim----output-R-tendance chimique de "qx" (kg/kg/s)
 c======================================================================
+      USE infotrac
+      use dimphy
+      IMPLICIT none
 #include "dimensions.h"
-#include "dimphy.h"
 #include "clesphys.h"
 #include "YOMCST.h"
Index: /trunk/libf/phytitan/printflag.F
===================================================================
--- /trunk/libf/phytitan/printflag.F	(revision 101)
+++ /trunk/libf/phytitan/printflag.F	(revision 102)
@@ -12,6 +12,6 @@
 
        REAL tabcntr0( 100 )
-       LOGICAL cycle_diurn0,soil_model0,new_oliq0,ok_orodr0
-       LOGICAL ok_orolf0,ok_limitvr0
+       LOGICAL cycle_diurn0,soil_model0,ok_orodr0
+       LOGICAL ok_orolf0,ok_gw_nonoro0
        LOGICAL ok_mensuel, ok_journe, ok_instan
        INTEGER radpas , radpas0
@@ -33,6 +33,6 @@
        PRINT 100
 
-c      PRINT 11, new_oliq, ok_orodr, ok_orolf   
-c      PRINT 100
+       PRINT 11, ok_orodr, ok_orolf, ok_gw_nonoro   
+       PRINT 100
 
        PRINT 12, nbapp_rad
@@ -55,15 +55,13 @@
         cycle_diurn0  = .FALSE.
         soil_model0   = .FALSE.
-        new_oliq0     = .FALSE.
         ok_orodr0     = .FALSE.
         ok_orolf0     = .FALSE.
-        ok_limitvr0   = .FALSE.
+        ok_gw_nonoro0 = .FALSE.
 
         IF( tabcntr0( 7 ).EQ. 1. )   cycle_diurn0 = .TRUE.
         IF( tabcntr0( 8 ).EQ. 1. )    soil_model0 = .TRUE.
-        IF( tabcntr0( 9 ).EQ. 1. )      new_oliq0 = .TRUE.
         IF( tabcntr0(10 ).EQ. 1. )      ok_orodr0 = .TRUE.
         IF( tabcntr0(11 ).EQ. 1. )      ok_orolf0 = .TRUE.
-        IF( tabcntr0(12 ).EQ. 1. )    ok_limitvr0 = .TRUE.
+        IF( tabcntr0(12 ).EQ. 1. )  ok_gw_nonoro0 = .TRUE.
 
         PRINT *,' $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
@@ -108,10 +106,4 @@
        ENDIF
 
-       IF( new_oliq0.AND..NOT.new_oliq.OR..NOT.new_oliq0.AND.
-     ,        new_oliq )     THEN
-        PRINT 16, new_oliq0, new_oliq
-        PRINT 100
-       ENDIF
-
        IF( ok_orodr0.AND..NOT.ok_orodr.OR..NOT.ok_orodr0.AND.
      ,        ok_orodr )     THEN
@@ -126,7 +118,7 @@
        ENDIF
 
-       IF( ok_limitvr0.AND..NOT.ok_limitvrai.OR..NOT.ok_limitvr0.
-     ,     AND.ok_limitvrai )     THEN
-        PRINT 18, ok_limitvr0, ok_limitvrai
+       IF( ok_gw_nonoro0.AND..NOT.ok_gw_nonoro.OR..NOT.ok_gw_nonoro0.
+     ,     AND.ok_gw_nonoro )     THEN
+        PRINT 18, ok_gw_nonoro0, ok_gw_nonoro
         PRINT 100
        ENDIF
@@ -140,6 +132,4 @@
      . , ok_instan = ',l3,8x,5(1H*) )
 
- 7     FORMAT(2x,5(1H*),15x,'      ok_limitvrai   = ',l3,16x,5(1h*) )
-
  8     FORMAT(2x,'*****             radpas    =                      ' ,
      , i4,6x,' *****')
@@ -149,6 +139,6 @@
 
 
- 11    FORMAT(2x,5(1H*),'  new_oliq = ',l3,3x,', Ok_orodr = ',
-     , l3,3x,', Ok_orolf = ',l3,3x,5(1H*) )
+ 11    FORMAT(2x,5(1H*),', Ok_orodr = ',l3,3x,', Ok_orolf = ',l3,3x,
+     ,    ' ok_gw_nonoro = ',l3,3x,5(1H*) )
 
 
@@ -165,11 +155,8 @@
      , /1x,10x,' startphy = ',l3,2x,' et  run.def = ',l3)
 
- 16    FORMAT(2x,'$$$$$$$$   Attention !!      new_oliq  different  sur',
-     , /1x,10x,' startphy = ',l3,2x,' et  run.def = ',l3)
-
  17    FORMAT(2x,'$$$$$$$$   Attention !!      ok_orolf  different  sur',
      , /1x,10x,' startphy = ',l3,2x,' et  run.def = ',l3)
 
- 18    FORMAT(2x,'$$$$$$$$   Attention !!  ok_limitvrai  different  sur',
+ 18    FORMAT(2x,'$$$$$$$$   Attention !!  ok_gw_nonoro  different  sur',
      , /1x,10x,' startphy = ',l3,2x,' et  run.def = ',l3)
 
Index: unk/libf/phytitan/raddim.h
===================================================================
--- /trunk/libf/phytitan/raddim.h	(revision 101)
+++ 	(revision )
@@ -1,5 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/raddim.48.32.h,v 1.1.1.1 2004/05/19 12:53:08 lmdzadmin Exp $
-!
-      INTEGER kdlon, kflev
-      PARAMETER (kdlon=149,kflev=klev)
Index: /trunk/libf/phytitan/radlwsw.F
===================================================================
--- /trunk/libf/phytitan/radlwsw.F	(revision 101)
+++ /trunk/libf/phytitan/radlwsw.F	(revision 102)
@@ -6,5 +6,4 @@
      .                  lwnet, swnet)
 c      
-      IMPLICIT none
 c======================================================================
 c Auteur(s): Z.X. Li (LMD/CNRS) date: 19960719
@@ -39,10 +38,10 @@
 
 c======================================================================
+      use dimphy
+      USE comgeomphy
+      IMPLICIT none
 #include "dimensions.h"
-#include "dimphy.h"
-#include "raddim.h"
 #include "YOMCST.h"
 #include "clesphys.h" 
-#include "comgeomphy.h"
 c
 c ARGUMENTS
Index: /trunk/libf/phytitan/radtitan.F
===================================================================
--- /trunk/libf/phytitan/radtitan.F	(revision 101)
+++ /trunk/libf/phytitan/radtitan.F	(revision 102)
@@ -1,6 +1,3 @@
        SUBROUTINE RADTITAN(p,nq,nmicro,ycomp)
-
-      IMPLICIT NONE
-
 
 c=======================================================================
@@ -34,12 +31,13 @@
 c   -------------
 
+      USE infotrac
+      use dimphy
+      USE comgeomphy
+      IMPLICIT NONE
 #include "dimensions.h"
-#include "dimphy.h"
 #include "clesphys.h"
 #include "microtab.h"
 #include "numchimrad.h"
-#include "comgeomphy.h"
 #include "YOMCST.h"
-#include "advtrac.h"     !! pour noms des traceurs
 
 c Pour le CRAY, les block data doivent etre declares external
@@ -297,15 +295,15 @@
 
           do i=1,nq
-             if (tnom(i).eq."CH4") then
+             if (tname(i).eq."CH4") then
                 iradch4=i
-             elseif (tnom(i).eq."C2H2") then
+             elseif (tname(i).eq."C2H2") then
                 iradc2h2=i
-             elseif (tnom(i).eq."C2H6") then
+             elseif (tname(i).eq."C2H6") then
                 iradc2h6=i
-             elseif (tnom(i).eq."HCN") then
+             elseif (tname(i).eq."HCN") then
                 iradhcn=i
-             elseif (tnom(i).eq."N2") then
+             elseif (tname(i).eq."N2") then
                 iradn2=i
-             elseif (tnom(i).eq."H2") then
+             elseif (tname(i).eq."H2") then
                 iradh2=i
              endif 
Index: /trunk/libf/phytitan/readstartphy.F
===================================================================
--- /trunk/libf/phytitan/readstartphy.F	(revision 101)
+++ /trunk/libf/phytitan/readstartphy.F	(revision 102)
@@ -9,11 +9,11 @@
      .           fder,radsol,
      .           tabcntr0)
-      IMPLICIT none
 c======================================================================
 c Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
 c Objet: Lecture de l'etat initial pour la physique
 c======================================================================
+      use dimphy
+      IMPLICIT none
 #include "dimensions.h"
-#include "dimphy.h"
 #include "netcdf.inc"
 #include "dimsoil.h"
Index: /trunk/libf/phytitan/sfluxv.F
===================================================================
--- /trunk/libf/phytitan/sfluxv.F	(revision 101)
+++ /trunk/libf/phytitan/sfluxv.F	(revision 102)
@@ -1,8 +1,7 @@
       SUBROUTINE SFLUXV(IPRINT,IG,dist_sol)
 
+      use dimphy
       IMPLICIT NONE
-
 #include "dimensions.h"
-#include "dimphy.h"
 #include "comorbit.h"
 
Index: /trunk/libf/phytitan/soil.F
===================================================================
--- /trunk/libf/phytitan/soil.F	(revision 101)
+++ /trunk/libf/phytitan/soil.F	(revision 102)
@@ -4,5 +4,4 @@
       SUBROUTINE soil(ptimestep, knon, ptsrf, ptsoil,
      s          pcapcal, pfluxgrd)
-      IMPLICIT NONE
 
 c=======================================================================
@@ -47,7 +46,8 @@
 c   -------------
 
+      use dimphy
+      IMPLICIT NONE
 #include "dimensions.h"
 #include "YOMCST.h"
-#include "dimphy.h"
 #include "dimsoil.h"
 #include "clesphys.h"
@@ -74,5 +74,5 @@
 c   ----------------------
       REAL dz1(nsoilmx),dz2(nsoilmx)
-      REAL zc(klon,nsoilmx),zd(klon,nsoilmx)
+      REAL,allocatable :: zc(:,:),zd(:,:)
       REAL lambda
       SAVE dz1,dz2,zc,zd,lambda
@@ -95,4 +95,6 @@
 
       IF (firstcall) THEN
+
+      allocate(zc(klon,nsoilmx),zd(klon,nsoilmx))
 
 c-----------------------------------------------------------------------
Index: /trunk/libf/phytitan/ustarhb.F
===================================================================
--- /trunk/libf/phytitan/ustarhb.F	(revision 101)
+++ /trunk/libf/phytitan/ustarhb.F	(revision 102)
@@ -3,5 +3,4 @@
 !
       SUBROUTINE ustarhb(knon,u,v,cd_m, ustar)
-      IMPLICIT none
 c======================================================================
 c Laurent Li (LMD/CNRS), le 30 septembre 1998
@@ -17,6 +16,7 @@
 c model. J. of Climate, vol. 6, 1825-1842.
 c======================================================================
+      use dimphy
+      IMPLICIT none
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 c
Index: /trunk/libf/phytitan/vdif_kcay.F
===================================================================
--- /trunk/libf/phytitan/vdif_kcay.F	(revision 101)
+++ /trunk/libf/phytitan/vdif_kcay.F	(revision 102)
@@ -5,8 +5,8 @@
      s   ,zlev,zlay,u,v,teta,cd,q2,q2diag,km,kn,ustar
      s   ,l_mix)
+c.......................................................................
+      use dimphy
       IMPLICIT NONE
-c.......................................................................
 #include "dimensions.h"
-#include "dimphy.h"
 c.......................................................................
 c
@@ -195,6 +195,4 @@
       PARAMETER (q2min=1.e-5)
       PARAMETER (q2max=1.E+2)
-      PARAMETER (nlay=klev)
-      PARAMETER (nlev=klev+1)
 c
       PARAMETER (
@@ -221,4 +219,7 @@
       save first
       data first/.true./
+
+      nlay=klev
+      nlev=klev+1
 c.......................................................................
 c  traitment des valeur de q2 en entree
Index: /trunk/libf/phytitan/writerestartphy.F
===================================================================
--- /trunk/libf/phytitan/writerestartphy.F	(revision 101)
+++ /trunk/libf/phytitan/writerestartphy.F	(revision 102)
@@ -6,11 +6,11 @@
      .           t_ancien)
 
-      IMPLICIT none
 
 c======================================================================
 c Ecriture de l'etat de redemarrage pour la physique (pour newstart)
 c======================================================================
+      use dimphy
+      IMPLICIT none
 #include "dimensions.h"
-#include "dimphy.h"
 #include "netcdf.inc"
 #include "dimsoil.h"
Index: /trunk/libf/phytitan/yamada.F
===================================================================
--- /trunk/libf/phytitan/yamada.F	(revision 101)
+++ /trunk/libf/phytitan/yamada.F	(revision 102)
@@ -5,8 +5,8 @@
      s   ,zlev,zlay,u,v,teta,cd,q2,km,kn,ustar
      s   ,l_mix)
+c.......................................................................
+      use dimphy
       IMPLICIT NONE
-c.......................................................................
 #include "dimensions.h"
-#include "dimphy.h"
 c.......................................................................
 c
@@ -46,13 +46,7 @@
       integer l_mix,ngrid
 
-
-      integer nlay,nlev
-      PARAMETER (nlay=klev)
-      PARAMETER (nlev=klev+1)
-
       logical first
       save first
       data first/.true./
-
 
       integer ig,k
@@ -98,6 +92,6 @@
 
 c  Correction d'un bug sauvage a verifier.
-c      do k=2,nlev
-      do k=2,nlay
+c      do k=2,klevp1
+      do k=2,klev
                                                           do ig=1,ngrid
          dz(ig,k)=zlay(ig,k)-zlay(ig,k-1)
Index: /trunk/libf/phytitan/yamada4.F
===================================================================
--- /trunk/libf/phytitan/yamada4.F	(revision 101)
+++ /trunk/libf/phytitan/yamada4.F	(revision 102)
@@ -6,8 +6,8 @@
      s   ,zlev,zlay,u,v,teta,cd,km,kn,kq,ustar
      s   ,iflag_pbl)
+c.......................................................................
+      use dimphy
       IMPLICIT NONE
-c.......................................................................
 #include "dimensions.h"
-#include "dimphy.h"
 c.......................................................................
 c
@@ -50,5 +50,5 @@
       REAL teta(klon,klev)
       REAL cd(klon)
-      REAL q2(klon,klev+1),qpre
+      REAL qpre
       REAL unsdz(klon,klev)
       REAL unsdzdec(klon,klev+1)
@@ -65,6 +65,4 @@
 
       integer nlay,nlev
-      PARAMETER (nlay=klev)
-      PARAMETER (nlev=klev+1)
 
       logical first
@@ -83,9 +81,9 @@
       real dtetadz(klon,klev+1)
       real m2cstat,mcstat,kmcstat
-      real l(klon,klev+1),l0(klon)
-      save l0
+      real l(klon,klev+1)
+      real,save,allocatable :: l0(:)
 c  ATTENTION! mis ici car j'ai enlevé q2 des arguments...
 c   sinon, c'est au-dessus que ça se passe...
-      save q2 
+      REAL,save,allocatable :: q2(:,:)
 
       real sq(klon),sqz(klon),zz(klon,klev+1)
@@ -99,8 +97,7 @@
       real fl,zzz,zl0,zq2,zn2
 
-      real rino(klon,klev+1),smyam(klon,klev),styam(klon,klev)
-     s  ,lyam(klon,klev),knyam(klon,klev)
-     s  ,w2yam(klon,klev),t2yam(klon,klev)
-      common/pbldiag/rino,smyam,styam,lyam,knyam,w2yam,t2yam
+c     real rino(klon,klev+1),smyam(klon,klev),styam(klon,klev)
+c    s  ,lyam(klon,klev),knyam(klon,klev)
+c    s  ,w2yam(klon,klev),t2yam(klon,klev)
 
       frif(ri)=0.6588*(ri+0.1776-sqrt(ri*ri-0.3221*ri+0.03156))
@@ -116,6 +113,13 @@
 
 c===================================
-c INITIALISATIONS (surtout pour k=1, à cause diagnostiques...)
+c INITIALISATIONS 
+      nlay=klev
+      nlev=klev+1
+
       if (first) then
+        allocate(l0(klon))
+        allocate(q2(klon,klevp1))
+
+c (surtout pour k=1, à cause diagnostiques...)
         q2 = 0.
       endif
@@ -462,22 +466,15 @@
 
 c     print*,'YAMADA4 1'
-c   Diagnostique pour stockage
-
-      rino=rif
-      smyam(:,1:klev)=sm(:,1:klev)
-      styam=sm(:,1:klev)*alpha(:,1:klev)
-      lyam(1:klon,1:klev)=l(:,1:klev)
-      knyam(1:klon,1:klev)=kn(:,1:klev)
 
 c   Estimations de w'2 et T'2 d'apres Abdela et McFarlane
 
-        if(1.eq.0)then
-      w2yam=q2(:,1:klev)*0.24
-     s    +lyam(:,1:klev)*5.17*kn(:,1:klev)*n2(:,1:klev)
-     s   /sqrt(q2(:,1:klev))
-
-      t2yam=9.1*kn(:,1:klev)*dtetadz(:,1:klev)**2/sqrt(q2(:,1:klev))
-     s  *lyam(:,1:klev)
-	endif
+c     if(1.eq.0)then
+c      w2yam=q2(:,1:klev)*0.24
+c    s    +lyam(:,1:klev)*5.17*kn(:,1:klev)*n2(:,1:klev)
+c    s   /sqrt(q2(:,1:klev))
+c
+c      t2yam=9.1*kn(:,1:klev)*dtetadz(:,1:klev)**2/sqrt(q2(:,1:klev))
+c    s  *lyam(:,1:klev)
+c     endif
 
 c     print*,'OKFIN'
Index: /trunk/libf/phytitan/zenang.F
===================================================================
--- /trunk/libf/phytitan/zenang.F	(revision 101)
+++ /trunk/libf/phytitan/zenang.F	(revision 102)
@@ -2,5 +2,4 @@
       SUBROUTINE zenang(longi,gmtime,pdtrad,lat,long,
      s                  pmu0,frac)
-      IMPLICIT none
 c=============================================================
 c Auteur : O. Boucher (LMD/CNRS)
@@ -25,6 +24,7 @@
 c frac-----OUTPUT: ensoleillement moyen entre gmtime et gmtime+(pdtrad/RDAY)
 c================================================================
+      use dimphy
+      IMPLICIT none
 #include "dimensions.h"
-#include "dimphy.h"
 #include "YOMCST.h"
 #include "comorbit.h"
Index: unk/libf/phyvenus/comgeomphy.h
===================================================================
--- /trunk/libf/phyvenus/comgeomphy.h	(revision 101)
+++ 	(revision )
@@ -1,9 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/comgeomphy.h,v 1.1.1.1 2004/05/19 12:53:09 lmdzadmin Exp $
-!
-c
-c
-c Common de passage de la geometrie de la dynamique a la physique
-      real airephy(klon),cuphy(klon),cvphy(klon)
-      REAL rlatd(klon), rlond(klon)
-      common/comgeomphy/airephy,cuphy,cvphy,rlatd, rlond
Index: unk/libf/phyvenus/dimphy.h
===================================================================
--- /trunk/libf/phyvenus/dimphy.h	(revision 101)
+++ 	(revision )
@@ -1,9 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/dimphy.h,v 1.1.1.1 2004/05/19 12:53:07 lmdzadmin Exp $
-!
-c-----------------------------------------------------------------------
-      INTEGER KIDIA, KFDIA, KLON, KLEV
-      PARAMETER (KIDIA=1,KFDIA=iim*(jjm-1)+2-1/jjm,
-     .           KLON=KFDIA-KIDIA+1,KLEV=llm)
-      REAL zmasq(KLON)
-      COMMON/terreoce/zmasq
Index: /trunk/libf/phyvenus/iniphysiq.F
===================================================================
--- /trunk/libf/phyvenus/iniphysiq.F	(revision 101)
+++ /trunk/libf/phyvenus/iniphysiq.F	(revision 102)
@@ -10,7 +10,4 @@
      $           prad,pg,pr,pcpp)
 
-      use dimphy
-      USE comgeomphy
-      IMPLICIT NONE
 c
 c=======================================================================
@@ -46,4 +43,7 @@
 c   -------------
  
+      use dimphy
+      USE comgeomphy
+      IMPLICIT NONE
 #include "dimensions.h"
 
Index: /trunk/libf/phyvenus/printflag.F
===================================================================
--- /trunk/libf/phyvenus/printflag.F	(revision 101)
+++ /trunk/libf/phyvenus/printflag.F	(revision 102)
@@ -119,6 +119,6 @@
 
 
- 11    FORMAT(2x,5(1H*),', Ok_orodr = ',
-     , l3,3x,', Ok_orolf = ',l3,3x,5(1H*) )
+ 11    FORMAT(2x,5(1H*),', Ok_orodr = ',l3,3x,', Ok_orolf = ',l3,3x,
+     ,    ' ok_gw_nonoro = ',l3,3x,5(1H*) )
 
 
Index: unk/libf/phyvenus/raddim.160.98.h
===================================================================
--- /trunk/libf/phyvenus/raddim.160.98.h	(revision 101)
+++ 	(revision )
@@ -1,5 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/raddim.160.98.h,v 1.1.1.1 2004/05/19 12:53:08 lmdzadmin Exp $
-!
-      INTEGER kdlon, kflev
-      PARAMETER (kdlon=78,kflev=klev) ! 78*199
Index: unk/libf/phyvenus/raddim.192.143.h
===================================================================
--- /trunk/libf/phyvenus/raddim.192.143.h	(revision 101)
+++ 	(revision )
@@ -1,5 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/raddim.192.143.h,v 1.1.1.1 2004/05/19 12:53:07 lmdzadmin Exp $
-!
-       INTEGER kdlon, kflev
-      PARAMETER (kdlon=10,kflev=klev)
Index: unk/libf/phyvenus/raddim.32.24.h
===================================================================
--- /trunk/libf/phyvenus/raddim.32.24.h	(revision 101)
+++ 	(revision )
@@ -1,5 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/raddim.32.24.h,v 1.1.1.1 2004/05/19 12:53:08 lmdzadmin Exp $
-!
-      INTEGER kdlon, kflev
-      PARAMETER (kdlon=klon,kflev=klev)
Index: unk/libf/phyvenus/raddim.48.32.h
===================================================================
--- /trunk/libf/phyvenus/raddim.48.32.h	(revision 101)
+++ 	(revision )
@@ -1,5 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/raddim.48.32.h,v 1.1.1.1 2004/05/19 12:53:08 lmdzadmin Exp $
-!
-      INTEGER kdlon, kflev
-      PARAMETER (kdlon=149,kflev=klev)
Index: unk/libf/phyvenus/raddim.72.46.h
===================================================================
--- /trunk/libf/phyvenus/raddim.72.46.h	(revision 101)
+++ 	(revision )
@@ -1,5 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/raddim.72.46.h,v 1.1.1.1 2004/05/19 12:53:08 lmdzadmin Exp $
-!
-      INTEGER kdlon, kflev
-      PARAMETER (kdlon=1621,kflev=klev)
Index: unk/libf/phyvenus/raddim.96.72.h
===================================================================
--- /trunk/libf/phyvenus/raddim.96.72.h	(revision 101)
+++ 	(revision )
@@ -1,5 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/raddim.96.72.h,v 1.1.1.1 2004/05/19 12:53:08 lmdzadmin Exp $
-!
-      INTEGER kdlon, kflev
-      PARAMETER (kdlon=487,kflev=klev)
Index: unk/libf/phyvenus/raddim.defaut.h
===================================================================
--- /trunk/libf/phyvenus/raddim.defaut.h	(revision 101)
+++ 	(revision )
@@ -1,5 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/raddim.defaut.h,v 1.1.1.1 2004/05/19 12:53:09 lmdzadmin Exp $
-!
-      INTEGER kdlon, kflev
-      PARAMETER (kdlon=klon,kflev=klev)
Index: unk/libf/phyvenus/raddim.h
===================================================================
--- /trunk/libf/phyvenus/raddim.h	(revision 101)
+++ 	(revision )
@@ -1,5 +1,0 @@
-!
-! $Header: /home/cvsroot/LMDZ4/libf/phylmd/raddim.48.32.h,v 1.1.1.1 2004/05/19 12:53:08 lmdzadmin Exp $
-!
-      INTEGER kdlon, kflev
-      PARAMETER (kdlon=149,kflev=klev)
