Index: /trunk/chantiers/commit_importants.log
===================================================================
--- /trunk/chantiers/commit_importants.log	(revision 50)
+++ /trunk/chantiers/commit_importants.log	(revision 51)
@@ -301,3 +301,31 @@
 - Avec ces modifs, la compilation marche sans physque, avec ou sans ioipsl
   et avec la physique terrestre phylmd.
- 
+
+********************************************
+**** commit_v51 (et quelques autres...) ****
+********************************************
+
+29/01/2011 --- A. SPIGA
+
+Nouveautes LMD_MM_MARS - LMD_LES_MARS - LMDZ.MARS [dossier mars]
+
+Ce commit fait le bilan des chantiers d'uniformisation des diverses sources du mesochelle
+[tous les commits par "spiga" depuis le debut sont relatifs a ces chantiers]
+
+Chantiers acheves avec compilation et execution testees:
+
+1. le meme coeur dynamique peut maintenant etre compile avec ancienne ou nouvelle physique; option -p sur makemeso.
+
+2. la meme interface module_lmd_driver est partagee par LMD_MM_MARS et LMD_LES_MARS
+
+3. le meme libf est partage entre modele mesoechelle et gcm
+   les quelques routines qui different sont notees meso_ et integrees aux sources du gcm
+
+=====> De 4 version du modele + 1 version de la nouvelle physique, on passe a une seule version du modele
+meso-echelle dont la physique est automatiquement indexee sur les modifications effectuees dans le GCM 
+[quelques routines mesoechelle mises a part]
+
+* Reste a tester LMD_LES_MARS avec la nouvelle physique. Modifications mineures attendues suite au point 2. 
+* Reste a tester le nesting qui devrait induire des modifications dues a des lignes trop longues... ce sera transparent du point de vue GCM.
+
+
Index: /trunk/mars/libf/phymars/callradite.F
===================================================================
--- /trunk/mars/libf/phymars/callradite.F	(revision 50)
+++ /trunk/mars/libf/phymars/callradite.F	(revision 51)
@@ -263,8 +263,8 @@
 c        PLEASE MAKE SURE that you set up the right number of
 c          scatterers in dimradmars.h (naerkind);
-c        name_iaer(1) = "dust_conrath"
-         name_iaer(1) = "dust_doubleq"
-c        name_iaer(2) = "dust_submicron"
-         name_iaer(2) = "h2o_ice"
+         name_iaer(1) = "dust_conrath"
+c        name_iaer(1) = "dust_doubleq"
+cc        name_iaer(2) = "dust_submicron"
+c         name_iaer(2) = "h2o_ice"
 c        ----------------------------------------------------------
 
Index: /trunk/mars/libf/phymars/dimradmars.h
===================================================================
--- /trunk/mars/libf/phymars/dimradmars.h	(revision 50)
+++ /trunk/mars/libf/phymars/dimradmars.h	(revision 51)
@@ -15,7 +15,7 @@
       INTEGER  NFLEV,NDLON,NDLO2,ndomainsz
 
-      parameter (ndomainsz=ngridmx)
-!     parameter (ndomainsz=(ngridmx-1)/20 + 1)
-!     parameter (ndomainsz=(ngridmx-1)/5 + 1) 
+!      parameter (ndomainsz=ngridmx)
+      parameter (ndomainsz=(ngridmx-1)/20 + 1)
+!      parameter (ndomainsz=(ngridmx-1)/5 + 1) 
 
       parameter (NFLEV=nlayermx,NDLON=ndomainsz) ! avec decoupage
@@ -27,5 +27,6 @@
 ! (ex: naerkind=2 if you use one dust mode and active ice ...)
       integer naerkind
-      parameter (naerkind=2)
+      parameter (naerkind=1)
+!      parameter (naerkind=2)
 
 ! Reference wavelengths used to compute reference optical depth (m)
@@ -76,5 +77,5 @@
 
       INTEGER, PARAMETER :: nsizemax = 60
-!     INTEGER, PARAMETER :: nsizemax = 1
+!      INTEGER, PARAMETER :: nsizemax = 1
 
 ! Various initialisation for LW radiative code
Index: /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/dyn_em/module_initialize_real.F
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/dyn_em/module_initialize_real.F	(revision 50)
+++ /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/dyn_em/module_initialize_real.F	(revision 51)
@@ -330,4 +330,7 @@
                                        !! grid%mars_isoil(i,k+1,j) is defined to large-scale value grid%em_isoil_gc   
                    ENDIF        
+                ENDIF
+                IF (grid%mars_tsoil(i,k,j) .lt. 20.) THEN  !!! une securite pour les anciens diagfi qui n'ont que 10 niveaux
+                   IF (k .ne. 1) grid%mars_tsoil(i,k,j) = grid%mars_tsoil(i,k-1,j)
                 ENDIF
                 !!!!!!!!!!!!!!!!! DONE in soil_setting.F
Index: /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/assim_aeropacity.F
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/assim_aeropacity.F	(revision 51)
+++ /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/assim_aeropacity.F	(revision 51)
@@ -0,0 +1,464 @@
+      SUBROUTINE aeropacity(ngrid,nlayer,nq,zday,pplay,pplev,ls,pq,
+     &    tauref,tau,aerosol,reffrad,
+     &    QREFvis3d,QREFir3d,omegaREFvis3d,omegaREFir3d)
+                                                   
+! to use  'getin'
+      USE ioipsl_getincom 
+       IMPLICIT NONE
+c=======================================================================
+c   subject:
+c   --------
+c   Computing aerosol optical depth in each gridbox.
+c
+c   author: F.Forget 
+c   ------
+c   update F. Montmessin (water ice scheme) 
+c      and S. Lebonnois (12/06/2003) compatibility dust/ice/chemistry
+c   update J.-B. Madeleine 2008-2009:
+c       - added 3D scattering by aerosols;
+c       - dustopacity transferred from physiq.F to callradite.F,
+c           and renamed into aeropacity.F;
+c   
+c   input:
+c   ----- 
+c   ngrid             Number of gridpoint of horizontal grid
+c   nlayer            Number of layer
+c   nq                Number of tracer
+c   zday                  Date (time since Ls=0, in martian days)
+c   ls                Solar longitude (Ls) , radian
+c   pplay,pplev       pressure (Pa) in the middle and boundary of each layer
+c   pq                Dust mixing ratio (used if tracer =T and active=T).
+c   reffrad(ngrid,nlayer,naerkind)  Aerosol effective radius
+c   QREFvis3d(ngridmx,nlayermx,naerkind) \ 3d extinction coefficients
+c   QREFir3d(ngridmx,nlayermx,naerkind)  / at reference wavelengths;
+c   omegaREFvis3d(ngridmx,nlayermx,naerkind) \ 3d single scat. albedo
+c   omegaREFir3d(ngridmx,nlayermx,naerkind)  / at reference wavelengths;
+c
+c   output:
+c   -------
+c   tauref       Prescribed mean column optical depth at 700 Pa 
+c   tau          Column total visible dust optical depth at each point
+c   aerosol      aerosol(ig,l,1) is the dust optical
+c                depth in layer l, grid point ig
+
+c
+c=======================================================================
+#include "dimensions.h"
+#include "dimphys.h"
+#include "callkeys.h"
+#include "comcstfi.h"
+#include "comgeomfi.h"
+#include "dimradmars.h"
+#include "yomaer.h"
+#include "tracer.h"
+#include "planete.h"
+
+c-----------------------------------------------------------------------
+c
+c    Declarations :
+c    --------------
+c
+c    Input/Output
+c    ------------
+      INTEGER ngrid,nlayer,nq
+
+      REAL ls,zday,expfactor    
+      REAL pplev(ngrid,nlayer+1),pplay(ngrid,nlayer)
+      REAL pq(ngrid,nlayer,nq)
+      REAL tauref(ngrid), tau(ngrid,naerkind)
+      REAL aerosol(ngrid,nlayer,naerkind)
+      REAL reffrad(ngrid,nlayer,naerkind)
+      REAL QREFvis3d(ngridmx,nlayermx,naerkind)
+      REAL QREFir3d(ngridmx,nlayermx,naerkind)
+      REAL omegaREFvis3d(ngridmx,nlayermx,naerkind)
+      REAL omegaREFir3d(ngridmx,nlayermx,naerkind)
+c
+c    Local variables :
+c    -----------------
+      INTEGER l,ig,iq,i,j
+      INTEGER iaer           ! Aerosol index
+      real topdust(ngridmx)
+      real zlsconst, zp
+      real taueq,tauS,tauN
+      real r0,reff,coefsize
+c     Mean Qext(vis)/Qext(ir) profile
+      real msolsir(nlayermx,naerkind)
+c     Mean Qext(ir)/Qabs(ir) profile
+      real mqextsqabs(nlayermx,naerkind)
+c     Variables used when multiple particle sizes are used
+c       for dust or water ice particles in the radiative transfer
+c       (see callradite.F for more information).
+      REAL taudusttmp(ngridmx)! Temporary dust opacity
+                               !   used before scaling
+      REAL taudustvis(ngridmx) ! Dust opacity after scaling
+      REAL taudusttes(ngridmx) ! Dust opacity at IR ref. wav. as
+                               !   "seen" by the GCM.
+      REAL taucloudvis(ngridmx)! Cloud opacity at visible
+                               !   reference wavelength
+      REAL taucloudtes(ngridmx)! Cloud opacity at infrared
+                               !   reference wavelength using
+                               !   Qabs instead of Qext
+                               !   (direct comparison with TES)
+c
+c   local saved variables
+c   ---------------------
+
+      REAL topdust0(ngridmx) 
+      SAVE topdust0
+
+      LOGICAL firstcall
+      DATA firstcall/.true./
+      SAVE firstcall
+
+! indexes of water ice and dust tracers:
+      INTEGER,SAVE :: nqdust(nqmx) ! to store the indexes of dust tracers
+      INTEGER,SAVE :: i_ice=0  ! water ice
+      CHARACTER(LEN=20) :: tracername ! to temporarly store text
+
+      call zerophys(ngrid*naerkind,tau)
+
+! identify tracers
+
+      IF (firstcall) THEN
+        ! identify tracers which are dust
+        i=0
+        DO iq=1,nq
+          tracername=noms(iq)
+          IF (tracername(1:4).eq."dust") THEN
+          i=i+1
+          nqdust(i)=iq
+          ENDIF
+        ENDDO
+
+        IF (water.AND.activice) THEN
+          i_ice=igcm_h2o_ice
+          write(*,*) "aeropacity: i_ice=",i_ice
+        ENDIF
+
+c       altitude of the top of the aerosol layer (km) at Ls=2.76rad:
+c       in the Viking year scenario
+        DO ig=1,ngrid
+            topdust0(ig)=60. -22.*SIN(lati(ig))**2
+        END DO
+
+c       typical profile of solsir and (1-w)^(-1):
+        call zerophys(nlayer*naerkind,msolsir)
+        call zerophys(nlayer*naerkind,mqextsqabs)
+        WRITE(*,*) "Typical profiles of solsir and Qext/Qabs(IR):"
+        DO iaer = 1, naerkind ! Loop on aerosol kind
+          WRITE(*,*) "Aerosol # ",iaer
+          DO l=1,nlayer
+            DO ig=1,ngridmx
+              msolsir(l,iaer)=msolsir(l,iaer)+
+     &              QREFvis3d(ig,l,iaer)/
+     &              QREFir3d(ig,l,iaer)
+              mqextsqabs(l,iaer)=mqextsqabs(l,iaer)+
+     &              (1.E0-omegaREFir3d(ig,l,iaer))**(-1)
+            ENDDO
+            msolsir(l,iaer)=msolsir(l,iaer)/REAL(ngridmx)
+            mqextsqabs(l,iaer)=mqextsqabs(l,iaer)/REAL(ngridmx)
+          ENDDO
+          WRITE(*,*) "solsir: ",msolsir(:,iaer)
+          WRITE(*,*) "Qext/Qabs(IR): ",mqextsqabs(:,iaer)
+        ENDDO
+
+!       load value of tauvis from callphys.def (if given there,
+!       otherwise default value read from starfi.nc file will be used)
+        call getin("tauvis",tauvis)
+
+        firstcall=.false.
+
+      END IF
+
+      DO iaer = 1, naerkind ! Loop on aerosol kind
+c     --------------------------------------------
+        aerkind: SELECT CASE (iaer)
+c==================================================================
+        CASE(1) aerkind                             ! Dust (iaer=1)
+c==================================================================
+
+c  -------------------------------------------------------------
+c     1) Prescribed dust  (if tracer=F or active=F)
+c  -------------------------------------------------------------
+      IF ((.not.tracer) .or. (.not.active)) THEN 
+
+c       Vertical column optical depth at 700.Pa 
+c       ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        IF(iaervar.eq.1) THEN 
+           do ig=1, ngridmx
+            tauref(ig)=max(tauvis,1.e-9) ! tauvis=cste (set in callphys.def
+                                         ! or read in starfi
+          end do
+        ELSE IF (iaervar.eq.2) THEN   ! << "Viking" Scenario>>
+
+          tauref(1) = 0.7+.3*cos(ls+80.*pi/180.) ! like seen by VL1
+          do ig=2,ngrid
+            tauref(ig) = tauref(1)
+          end do
+
+        ELSE IF (iaervar.eq.3) THEN  ! << "MGS" scenario >>
+
+           taueq= 0.2 +(0.5-0.2) *(cos(0.5*(ls-4.363)))**14
+           tauS= 0.1 +(0.5-0.1)  *(cos(0.5*(ls-4.363)))**14
+           tauN = 0.1
+c	   if (peri_day.eq.150) then
+c	     tauS=0.1
+c	     tauN=0.1 +(0.5-0.1)  *(cos(0.5*(ls+pi-4.363)))**14
+c	     taueq= 0.2 +(0.5-0.2) *(cos(0.5*(ls+pi-4.363)))**14
+c           endif
+           do ig=1,ngrid/2  ! Northern hemisphere
+             tauref(ig)= tauN +
+     &       (taueq-tauN)*0.5*(1+tanh((45-lati(ig)*180./pi)*6/60))
+           end do
+           do ig=ngrid/2+1, ngridmx  ! Southern hemisphere
+             tauref(ig)= tauS +
+     &       (taueq-tauS)*0.5*(1+tanh((45+lati(ig)*180./pi)*6/60))
+           end do
+
+        ELSE IF ((iaervar.eq.4).or.
+     &           ((iaervar.ge.24).and.(iaervar.le.26)))
+     &       THEN  ! << "TES assimilated dust scenarios >>
+           call readtesassim(ngrid,nlayer,zday,pplev,tauref)
+
+        ELSE IF (iaervar.eq.5) THEN   ! << Escalier Scenario>>
+c         tauref(1) = 0.2
+c         if ((ls.ge.210.*pi/180.).and.(ls.le.330.*pi/180.))
+c    &                              tauref(1) = 2.5
+          tauref(1) = 2.5
+          if ((ls.ge.30.*pi/180.).and.(ls.le.150.*pi/180.))
+     &                              tauref(1) = .2
+          do ig=2,ngrid
+            tauref(ig) = tauref(1)
+          end do
+
+        ELSE IF (iaervar.gt.99) THEN  ! << input netcdf file >>
+c*************WRF
+c
+c 2. customized dust opacity field
+c   ex: from assimilation
+       call meso_readtesassim(ngrid,nlayer,zday,pplev,tauref,
+     . iaervar)
+c
+c*************WRF
+
+        ELSE
+          stop 'problem with iaervar in aeropacity.F'
+        ENDIF
+
+c       Altitude of the top of the dust layer
+c       ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        zlsconst=SIN(ls-2.76)
+        if (iddist.eq.1) then
+          do ig=1,ngrid
+             topdust(ig)=topdustref         ! constant dust layer top
+          end do
+
+        else if (iddist.eq.2) then          ! "Viking" scenario
+          do ig=1,ngrid
+            topdust(ig)=topdust0(ig)+18.*zlsconst
+          end do
+
+        else if(iddist.eq.3) then         !"MGS" scenario
+          do ig=1,ngrid
+            topdust(ig)=60.+18.*zlsconst
+     &                -(32+18*zlsconst)*sin(lati(ig))**4
+     &                 - 8*zlsconst*(sin(lati(ig)))**5
+          end do
+        endif
+
+
+c       Optical depth in each layer :
+c       ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        if(iddist.ge.1) then
+          expfactor=0.
+          CALL zerophys(ngrid,taudusttmp)
+          DO l=1,nlayer
+            DO ig=1,ngrid
+c             Typical mixing ratio profile 
+              if(pplay(ig,l).gt.700.
+     $                        /(988.**(topdust(ig)/70.))) then
+                zp=(700./pplay(ig,l))**(70./topdust(ig))
+                 expfactor=max(exp(0.007*(1.-max(zp,1.))),1.e-3)
+              else    
+                expfactor=1.e-3
+              endif
+c             Vertical scaling function
+              aerosol(ig,l,iaer)= (pplev(ig,l)-pplev(ig,l+1)) *
+     &          expfactor *
+     &          QREFvis3d(ig,l,iaer) / QREFvis3d(ig,1,iaer)
+c             Scaling factor
+              taudusttmp(ig)=taudusttmp(ig)+aerosol(ig,l,iaer)
+            ENDDO
+          ENDDO
+
+c         Rescaling each layer to reproduce the choosen (or
+c           assimilated) dust extinction opacity at visible
+c           reference wavelength, which is originally scaled
+c           to an equivalent 700Pa pressure surface.
+          DO l=1,nlayer
+            DO ig=1,ngrid
+              aerosol(ig,l,iaer) = tauref(ig) *
+     &                     pplev(ig,1) / 700.E0 *
+     &                     aerosol(ig,l,iaer) / taudusttmp(ig)
+            ENDDO
+          ENDDO
+
+          CALL zerophys(ngrid,taudustvis)
+          CALL zerophys(ngrid,taudusttes)
+          DO l=1,nlayer
+            DO ig=1,ngrid
+              taudustvis(ig) = taudustvis(ig) + aerosol(ig,l,iaer)
+              taudusttes(ig) = taudusttes(ig) + aerosol(ig,l,iaer)*
+     &          QREFir3d(ig,l,iaer) / QREFvis3d(ig,l,iaer)*
+     &          ( 1.E0 - omegaREFir3d(ig,l,iaer) )
+            ENDDO
+          ENDDO
+
+c         Outputs
+          IF (ngrid.NE.1) THEN
+            CALL WRITEDIAGFI(ngridmx,'taudustTES','dust abs IR',
+     &        ' ',2,taudusttes)
+            CALL wstats(ngridmx,'taudustTES','dust abs IR',
+     &        ' ',2,taudusttes)
+          ELSE
+            CALL writeg1d(ngrid,1,taudusttes,'taudusttes','NU')
+          ENDIF
+
+c     changement dans le calcul de la distribution verticale          
+c     dans le cas des scenarios de poussieres assimiles
+c        if (iaervar.eq.4) THEN  ! TES
+c        call zerophys(ngrid*naerkind,tau)
+c
+c        do l=1,nlayer
+c           do ig=1,ngrid
+c                tau(ig,1)=tau(ig,1)+ aerosol(ig,l,1)
+c           end do
+c        end do
+c        do l=1,nlayer
+c           do ig=1,ngrid
+c               aerosol(ig,l,1)=aerosol(ig,l,1)*tauref(ig)/tau(ig,1)
+c     $     *(pplev(ig,1)/700)
+c           end do
+c        end do
+c        endif
+cccccccccccccccccccccccccccccccccccccccccccccccccccccccccc       
+        else if(iddist.eq.0) then   
+c         old dust vertical distribution function (pollack90)
+          DO l=1,nlayer
+             DO ig=1,ngrid
+                zp=700./pplay(ig,l)
+                aerosol(ig,l,1)= tauref(ig)/700. *
+     s           (pplev(ig,l)-pplev(ig,l+1))
+     s           *max( exp(.03*(1.-max(zp,1.))) , 1.E-3 )
+             ENDDO
+          ENDDO
+        end if
+
+c  ---------------------------------------------------------------------
+c     2) Transported radiatively active dust  (if tracer=T and active=T)
+c ----------------------------------------------------------------------
+      ELSE  IF ((tracer) .and. (active)) THEN 
+c     The dust opacity is computed from q
+
+c     a) "doubleq" technique (transport of mass and number mixing ratio)
+c        ~~~~~~~~~~~~~~~~~~~
+       if(doubleq) then 
+
+         call zerophys(ngrid*nlayer*naerkind,aerosol)
+
+c        Computing effective radius :
+         do  l=1,nlayer
+           do ig=1, ngrid
+              r0=
+     &        (r3n_q*pq(ig,l,1)/max(pq(ig,l,2),0.01))**(1./3.)
+              r0=min(max(r0,1.e-10),500.e-6)
+              reff= ref_r0 * r0
+cc           If  reff is small, the transported dust mean Qext
+c            is reduced from the reference dust Qext by a factor "coefsize"
+
+             coefsize=min(max(2.52e6*reff-0.043 ,0.)   ,1.)
+
+cc              It is added 1.e-8 to pq to avoid low
+
+                aerosol(ig,l,1)=aerosol(ig,l,1)+ 1.E-8 +
+     &         ( 0.75*Qext(1)*coefsize/(rho_dust*reff))
+     &          * (pq(ig,l,nqdust(1)))*
+c               only one dust bin to use with doubleq
+     &          (pplev(ig,l)-pplev(ig,l+1))/g
+           end do
+         end do
+         call zerophys(ngrid,tauref)
+
+c     b) Size bin technique (each aerosol can contribute to opacity))
+c        ~~~~~~~~~~~~~~~~~~
+       else
+c        The dust opacity is computed from q
+         call zerophys(ngrid*nlayer*naerkind,aerosol)
+         do iq=1,dustbin
+           do l=1,nlayer
+              do ig=1,ngrid
+cc               qextrhor(iq) is  (3/4)*Qext/(rho*reff)
+cc               It is added 1.e-8 to pq to avoid low
+                 aerosol(ig,l,1)=aerosol(ig,l,1)+
+     &           qextrhor(nqdust(iq))*(pq(ig,l,nqdust(iq))+1.e-8)*
+     &           (pplev(ig,l)-pplev(ig,l+1))/g
+              end do
+           end do
+         end do
+         call zerophys(ngrid,tauref)
+       end if  ! (doubleq)
+      END IF   ! (dust scenario)
+
+
+c==================================================================
+        CASE(2) aerkind               ! Water ice crystals (iaer=2)
+c==================================================================
+
+      IF (water.AND.activice) THEN
+c       1. Initialization
+        CALL zerophys(ngrid*nlayer,aerosol(1,1,iaer))
+        CALL zerophys(ngrid,taucloudvis)
+        CALL zerophys(ngrid,taucloudtes)
+c       2. Opacity calculation
+        DO ig=1, ngrid
+          DO l=1,nlayer
+            aerosol(ig,l,iaer) =
+     &        (  0.75 * QREFvis3d(ig,l,iaer) /
+     &        ( rho_ice * reffrad(ig,l,iaer) )  ) *
+     &        ( pq(ig,l,i_ice) + 1.E-8 ) *
+     &        ( pplev(ig,l) - pplev(ig,l+1) ) / g
+            taucloudvis(ig) = taucloudvis(ig) + aerosol(ig,l,iaer)
+            taucloudtes(ig) = taucloudtes(ig) + aerosol(ig,l,iaer)*
+     &        QREFir3d(ig,l,iaer) / QREFvis3d(ig,l,iaer) *
+     &        ( 1.E0 - omegaREFir3d(ig,l,iaer) )
+          ENDDO
+        ENDDO
+c       3. Outputs
+        IF (ngrid.NE.1) THEN
+          CALL WRITEDIAGFI(ngridmx,'tauTES','tauabs IR refwvl',
+     &      ' ',2,taucloudtes)
+          CALL wstats(ngridmx,'tauTES','tauabs IR refwvl',
+     &      ' ',2,taucloudtes)
+        ELSE
+          CALL writeg1d(ngrid,1,taucloudtes,'tautes','NU')
+        ENDIF
+      ENDIF
+
+c==================================================================
+        END SELECT aerkind
+
+c -----------------------------------------------------------------
+c Column integrated visible optical depth in each point
+c -----------------------------------------------------------------
+
+      do l=1,nlayer
+         do ig=1,ngrid
+               tau(ig,iaer) = tau(ig,iaer) + aerosol(ig,l,iaer)
+         end do
+      end do
+
+c     -----------------------------------
+      ENDDO ! iaer (loop on aerosol kind)
+
+      return
+      end 
Index: /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/assim_readtesassim.F90
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/assim_readtesassim.F90	(revision 51)
+++ /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/assim_readtesassim.F90	(revision 51)
@@ -0,0 +1,306 @@
+subroutine meso_readtesassim(ngrid,nlayer,zday,pplev,tauref,day_translate)
+
+! Reading of the dust assimilation file
+
+implicit none
+
+include "dimensions.h"
+include "dimphys.h"
+include "comgeomfi.h"
+include "netcdf.inc"
+include "datafile.h"
+
+integer, intent(in) :: ngrid,nlayer
+real, intent(in) :: zday ! date in martian days
+real, dimension(ngrid,nlayer+1), intent(in) :: pplev
+real, dimension(ngrid), intent(out) :: tauref
+
+! Local variables
+
+real :: realday
+integer nid,nvarid,ierr
+integer ltloop,lsloop,iloop,jloop,varloop,ig
+real, dimension(2) :: taubuf
+real tau1(4),tau
+real alt(4)
+integer latp(4),lonp(4)
+integer yinf,ysup,xinf,xsup,tinf,tsup
+real latinf,latsup,loninf,lonsup
+real latintmp,lonintmp,latdeg,londeg
+real colat,dlat,dlon,colattmp
+logical, save :: firstcall=.true.
+logical :: timeflag
+real radeg,pi
+integer :: timedim,londim,latdim
+real, dimension(:), allocatable, save :: lat,lon,time
+real, dimension(:,:,:), allocatable, save :: tautes
+integer, save :: timelen,lonlen,latlen
+INTEGER, external:: lnblnk
+
+!!WRF added
+INTEGER :: day_translate 
+
+pi=acos(-1.)
+radeg=180/pi
+realday=mod(zday,669.)
+
+if (firstcall) then
+   firstcall=.false.
+
+   !! Note: datafile() is defined in "datafile.h"
+   !print *,datafile(1:lnblnk(datafile))
+   !ierr=NF_OPEN (datafile(1:lnblnk(datafile))//"/input_totoptdep.nc",NF_NOWRITE,nid)
+
+!!-------------------------------------
+!! iaervar is used here to communicate the day starting the diagfi
+!! -- in addition, inpu_totoptdep.nc must be prepared
+day_translate=day_translate-1000
+ierr=NF_OPEN ("./input_totoptdep.nc",NF_NOWRITE,nid)
+!!-------------------------------------
+
+
+   IF (ierr.NE.NF_NOERR) THEN
+     !write(*,*)'Problem opening dust.nc (in phymars/readtesassim.F90)'
+     !write(*,*)'It should be in :',datafile(1:lnblnk(datafile)),'/'
+     !write(*,*)'1) You can change this directory address in '
+     !write(*,*)'   file phymars/datafile.h'
+     !write(*,*)'2) If necessary, dust.nc (and other datafiles)'
+     !write(*,*)'   can be obtained online on:'
+     !write(*,*)'   http://www.lmd.jussieu.fr/~forget/datagcm/datafile'
+     write(*,*)'Problem opening input_totoptdep.nc ...'
+     CALL ABORT
+   ENDIF
+
+   ierr=NF_INQ_DIMID(nid,"time",timedim)
+   ierr=NF_INQ_DIMLEN(nid,timedim,timelen)
+   ierr=NF_INQ_DIMID(nid,"lat",latdim)
+   ierr=NF_INQ_DIMLEN(nid,latdim,latlen)
+   ierr=NF_INQ_DIMID(nid,"lon",londim)
+   ierr=NF_INQ_DIMLEN(nid,londim,lonlen)
+
+   allocate(tautes(lonlen,latlen,timelen))
+   allocate(lat(latlen), lon(lonlen), time(timelen))
+
+   ierr = NF_INQ_VARID (nid, "totoptdep",nvarid)
+#ifdef NC_DOUBLE
+   ierr = NF_GET_VAR_DOUBLE(nid, nvarid, tautes)
+#else
+   ierr = NF_GET_VAR_REAL(nid, nvarid, tautes)
+#endif
+   IF (ierr .NE. NF_NOERR) THEN
+      PRINT*, "Error: Readtesassim <totoptdep> not found"
+      stop
+   ENDIF
+
+!!****WRF
+!! diagfi files feature the dust opacity divided by 700 Pa
+        tautes=tautes*700.
+!!****WRF
+
+
+   ierr = NF_INQ_VARID (nid, "time",nvarid)
+#ifdef NC_DOUBLE
+   ierr = NF_GET_VAR_DOUBLE(nid, nvarid, time)
+#else
+   ierr = NF_GET_VAR_REAL(nid, nvarid, time)
+#endif
+   IF (ierr .NE. NF_NOERR) THEN
+      PRINT*, "Error: Readtesassim <Time> not found"
+      stop
+   ENDIF
+
+!!****WRF
+PRINT *, "**** read assimilation dust"
+PRINT *, "prescribed starting day ...", day_translate
+time=time+day_translate
+!!****WRF
+
+   ierr = NF_INQ_VARID (nid, "lat",nvarid)
+#ifdef NC_DOUBLE
+   ierr = NF_GET_VAR_DOUBLE(nid, nvarid, lat)
+#else
+   ierr = NF_GET_VAR_REAL(nid, nvarid, lat)
+#endif
+   IF (ierr .NE. NF_NOERR) THEN
+      PRINT*, "Error: Readtesassim <latitude> not found"
+      stop
+   ENDIF
+
+   ierr = NF_INQ_VARID (nid, "lon",nvarid)
+#ifdef NC_DOUBLE
+   ierr = NF_GET_VAR_DOUBLE(nid, nvarid, lon)
+#else
+   ierr = NF_GET_VAR_REAL(nid, nvarid, lon)
+#endif
+   IF (ierr .NE. NF_NOERR) THEN
+      PRINT*, "Error: Readtesassim <longitude> not found"
+      stop
+   ENDIF
+
+endif ! of if (firstcall)
+
+do ig=1,ngrid
+
+! Find the four nearest points, arranged as follows:
+!                               1 2
+!                               3 4
+
+   latdeg=lati(ig)*radeg  ! latitude, in degrees
+   londeg=long(ig)*radeg  ! longitude, in degrees east
+   colat=90-latdeg        ! colatitude, in degrees
+
+ ! Find encompassing latitudes
+   if (colat<(90-lat(1))) then ! between north pole and lat(1)
+      ysup=1
+      yinf=1
+   else if (colat>=90-(lat(latlen))) then ! between south pole and lat(laten)
+      ysup=latlen
+      yinf=latlen
+   else ! general case
+      do iloop=2,latlen
+         if(colat<(90-lat(iloop))) then
+           ysup=iloop-1
+           yinf=iloop
+           exit
+         endif
+      enddo
+   endif
+
+   latinf=lat(yinf)
+   latsup=lat(ysup)
+
+
+ ! Find encompassing longitudes
+   ! Note: in input file, lon(1)=-180.
+   if (londeg>lon(lonlen)) then
+      xsup=1
+      xinf=lonlen
+      loninf=lon(xsup)
+      lonsup=180.0 ! since lon(1)=-180.0
+   else
+      do iloop=2,lonlen
+         if(londeg<lon(iloop)) then
+              xsup=iloop
+              xinf=iloop-1
+              exit
+         endif
+      enddo
+      loninf=lon(xinf)
+      lonsup=lon(xsup)
+   endif
+
+   if ((xsup.gt.lonlen).OR.(yinf.gt.latlen).OR.(xinf.lt.1)&
+     .OR.(ysup.lt.1)) then
+      write (*,*) "Readtesassim: SYSTEM ERROR on x or y in ts_gcm"
+      write (*,*) "xinf: ",xinf
+      write (*,*) "xsup: ",xsup
+      write (*,*) "yinf: ",yinf
+      write (*,*) "ysup: ",ysup
+      stop
+   endif
+
+!   loninf=lon(xinf)
+!   lonsup=lon(xsup)
+!   latinf=lat(yinf)
+!   latsup=lat(ysup)
+
+! The four neighbouring points are arranged as follows:
+!                               1 2
+!                               3 4
+
+   latp(1)=ysup
+   latp(2)=ysup
+   latp(3)=yinf
+   latp(4)=yinf
+
+   lonp(1)=xinf
+   lonp(2)=xsup
+   lonp(3)=xinf
+   lonp(4)=xsup
+
+! Linear interpolation on time, for all four neighbouring points
+
+   if ((realday<time(1)).or.(realday>time(timelen))) then
+      tinf=timelen
+      tsup=1
+      timeflag=.true.
+   else 
+      do iloop=2,timelen
+         if (realday<time(iloop)) then
+            tinf=iloop-1
+            tsup=iloop
+            exit
+         endif
+      enddo
+   endif
+
+! Bilinear interpolation on the four nearest points
+
+   if ((colat<(90-lat(1))).OR.(colat>(90-lat(latlen))).OR.(latsup==latinf)) then
+      dlat=0
+   else
+      dlat=((90-latinf)-colat)/((90-latinf)-(90-latsup))
+   endif
+
+   if (lonsup==loninf) then
+      dlon=0
+   else
+      dlon=(londeg-loninf)/(lonsup-loninf)
+   endif
+
+   do iloop=1,4
+      taubuf(1)=tautes(lonp(iloop),latp(iloop),tinf)
+      taubuf(2)=tautes(lonp(iloop),latp(iloop),tsup)
+      if (timeflag) then
+         if (realday>time(timelen)) then
+            tau1(iloop)=taubuf(1)+(taubuf(2)-taubuf(1))*(realday-time(tinf))/(time(tsup)+(669-time(tinf))) 
+         else
+            tau1(iloop)=taubuf(1)+(taubuf(2)-taubuf(1))*realday/(time(tsup)+(669-time(tinf)))
+         endif
+      else
+         tau1(iloop)=taubuf(1)+(taubuf(2)-taubuf(1))*(realday-time(tinf))/(time(tsup)-time(tinf))
+      endif
+      if (tau1(iloop)<0) then
+          write (*,*) "Readtesassim: SYSTEM ERROR on tau"
+          write (*,*) "utime ",realday
+          write (*,*) "time(tinf) ",time(tinf)
+          write (*,*) "time(tsup) ",time(tsup)
+          write (*,*) "tau1 ",taubuf(1)
+          write (*,*) "tau2 ",taubuf(2)
+          write (*,*) "tau ",tau1(iloop)
+          stop
+      endif
+   enddo
+   timeflag=.false.
+
+   if ((dlat>1).OR.(dlon>1) .OR. (dlat<0) .OR. (dlon<0)) then
+      write (*,*) "Readtesassim: SYSTEM ERROR on dlat or dlon in ts_gcm"
+      write (*,*) "dlat: ",dlat
+      write (*,*) "lat: ",latdeg
+      write (*,*) "dlon: ",dlon
+      write (*,*) "lon: ",londeg
+      stop
+   endif
+
+   tau= dlat*(dlon*(tau1(2)+tau1(3)-tau1(1)-tau1(4))+tau1(1)-tau1(3)) +dlon*(tau1(4)-tau1(3))+tau1(3)
+
+!!   tauref(ig)=tau*700/pplev(ig,1)*0.8
+    tauref(ig)=tau*0.825
+
+!-------------
+! from Oxford
+!-------------
+! + tuning ... TES is little bit too hot compared to radio-occultations 
+! 0.825 is the MY24 value
+!tauref(ig)=tau*700.*0.825
+!tauref(ig)=tau*700.
+!tauref(ig)=tau*700.*1.2
+!tauref(ig)=tau*700.*0.6
+!-------------------------------------------------
+! NB: multiply by 700. is now done elsewhere
+!-------------------------------------------------
+
+
+enddo
+
+end
Index: /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/diff_w_phys.sh
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/diff_w_phys.sh	(revision 51)
+++ /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/diff_w_phys.sh	(revision 51)
@@ -0,0 +1,21 @@
+#! /bin/bash
+
+
+logfile="diff_meso_gcm.log"
+wheregcm=$1
+
+\rm $logfile
+touch $logfile
+#echo '***** meso_inifis.F' >> $logfile
+diff ./meso_inifis.F      $1/inifis.F          >> $logfile
+#echo '***** meso_callkeys.h' >> $logfile 
+diff ./meso_callkeys.h    $1/callkeys.h        >> $logfile
+#echo '***** meso_dustlift.F' >> $logfile
+#diff ./meso_dustlift.F    $1/dustlift.F   >> $logfile
+#echo '***** meso_inifis.F' >> $logfile
+diff ./meso_inifis.F      $1/inifis.F          >> $logfile
+#echo '***** meso_newcondens.F' >> $logfile
+diff ./meso_newcondens.F  $1/newcondens.F      >> $logfile
+#echo '***** meso_physiq.F' >> $logfile
+diff ./meso_physiq.F      $1/physiq.F          >> $logfile
+diff ./meso_testphys1d.F  $1/testphys1d.F      >> $logfile
Index: /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/newgcm
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/newgcm	(revision 51)
+++ /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/newgcm	(revision 51)
@@ -0,0 +1,1 @@
+link /u/emlmd/LMDZ.MARS.BETA/libf/phymars
Index: /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/oldgcm
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/oldgcm	(revision 51)
+++ /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/oldgcm	(revision 51)
@@ -0,0 +1,1 @@
+link /donnees/aslmd/MODELES/WORKMODELS/LMDZ.MARS/libf/phymars
Index: /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/oldmeso
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/oldmeso	(revision 51)
+++ /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/oldmeso	(revision 51)
@@ -0,0 +1,1 @@
+link /donnees/aslmd/SVN/trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/libf/phymars
Index: /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/param_slope_full.F90
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/param_slope_full.F90	(revision 51)
+++ /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/param_slope_full.F90	(revision 51)
@@ -0,0 +1,248 @@
+SUBROUTINE param_slope_full( &
+  !
+  ! INPUTS
+  !
+               &  ls, localtime, latitude, taudust, albedo  &
+               &  ,theta_s, psi_s &
+               &  ,ftot_0  &
+  !
+  ! OUTPUTS
+  !
+               &  ,ftot &
+ )
+
+
+!!*****************************************************************************************
+!
+! SUBROUTINE:
+! param_slope
+! 
+!
+! PURPOSE: 
+! computes total solar irradiance on a given Martian slope 
+!
+!
+! INPUTS:
+! ls	 	aerocentric longitude (deg)
+! localtime	local true solar time (Martian hours)
+! latitude  	latitude (deg)
+! taudust	dust optical depth at reference wavelength 0.67 mic. 
+! albedo	spectrally integrated surface Lambertian reflection albedo 
+! theta_s	slope inclination angle (deg)
+!			0 is horizontal, 90 is vertical
+! phi_s		slope azimuth (deg) 
+!       		0 >> Northward
+!       		90 >> Eastward
+!       		180 >> Southward
+!       		270 >> Westward
+! ftot_0	spectrally integrated total irradiance on an horizontal surface (W/m2)
+!
+!
+! OUTPUTS:
+! ftot		spectrally integrated total irradiance on the slope (W/m2)
+!
+! REFERENCE: 
+! "Fast and accurate estimation of irradiance on Martian slopes"
+! A. Spiga & F. Forget
+! .....
+!
+! AUTHOR: 
+! A. Spiga (spiga@lmd.jussieu.fr)
+! March 2008
+!
+!!*****************************************************************************************
+
+IMPLICIT NONE
+
+!!
+!! INPUT
+!!
+REAL, INTENT(IN) :: ls, localtime, latitude, taudust, theta_s, psi_s, albedo, ftot_0  
+
+!!
+!! LOCAL
+!!
+REAL :: pi, deg2rad, dist_sol, cste_mars
+REAL, PARAMETER :: p = 1.510404          ! Semi-latus rectum of Martian elliptic orbit (AU)
+REAL, PARAMETER :: e = 9.3357898E-02     ! Eccentricity of Martian elliptic orbit
+REAL, PARAMETER :: t = 1.908231          ! Angle from Ls=0 to the perihelion (radian)
+REAL, PARAMETER :: so = 0.4256214        ! sin(Obliquity of Martian axis) 
+REAL :: rho, sdec, dec, cdec, csza, sza, ssza, psi0 
+REAL :: px, py
+REAL :: a
+REAL :: mu_s, sigma_s
+REAL :: fdir, fdir_0, fscat, fscat_0, fref
+REAL, DIMENSION(4,2) :: mat_M, mat_N, mat_T
+REAL, DIMENSION(2) :: g_vector
+REAL, DIMENSION(4) :: s_vector
+REAL :: ratio
+
+!!
+!! OUTPUT
+!!
+REAL, INTENT(OUT) :: ftot
+
+!!*****************************************************************************************
+
+!
+! Prerequisite
+!
+  pi = 2.*asin(1.)
+  deg2rad = pi/180.
+  if ((theta_s > 90.) .or. (theta_s < 0.)) then
+	print *, 'please set theta_s between 0 and 90', theta_s
+	stop 
+  endif
+
+!
+! Sun right ascension (radian)
+!
+  rho = pi*(1.0-localtime/12.0)
+
+!
+! Distance to sun (AU)
+!
+  dist_sol = p/(1.0+e*cos(deg2rad*Ls+t))   !! ellipse polar equation 
+ 
+!
+! Incident flux @ top of the atmosphere (Mars solar constant, W m-2)
+!
+  cste_mars=1370./(dist_sol*dist_sol)      !! 1370 W.m-2 is the solar constant at 1 AU.
+
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!
+!!! pour comparer avec spectres ESA
+!!!!!!!!!!!!!!!!!!!!!!!!!!!
+!cste_mars=cste_mars*0.92
+
+
+!
+! Sun declination (radian) [= subsolar point latitude]
+!
+  sdec = sin(deg2rad*Ls)*so
+  dec = asin(sdec)  
+  cdec = cos(dec)
+
+!
+! Solar Zenith angle (radian) 
+!
+  csza = sin(deg2rad*latitude)*sdec + cos(deg2rad*latitude)*cdec*cos(rho)
+  sza = acos(csza)
+  ssza = sin(sza)
+  if (csza < 0.01) then
+      !print *, 'sun below horizon' 
+      fdir_0=0.
+      fdir=0.
+      fscat_0=0.
+      fscat=0.   
+      fref=0.
+  else      
+
+!
+! 'Slope vs Sun' azimuth (radian)
+!
+  if ( ( (cdec*sin(rho)) .eq. 0.0 ) .and. ( ( sin(deg2rad*latitude)*cdec*cos(rho)-cos(deg2rad*latitude)*sdec ) .eq. 0.0 ) ) then
+    a = deg2rad*psi_s ! some compilator need specfying value for atan2(0,0)  
+  else
+    a = deg2rad*psi_s + atan2(cdec*sin(rho),sin(deg2rad*latitude)*cdec*cos(rho)-cos(deg2rad*latitude)*sdec)
+  end if
+  
+!
+! Cosine of slope-sun phase angle 
+!
+  mu_s = csza*cos(deg2rad*theta_s) - cos(a)*sin(deg2rad*theta_s)*sqrt(1-csza**2)
+  if (mu_s .le. 0.) mu_s=0.
+
+!
+! Sky-view factor
+!
+  sigma_s=0.5*(1.+cos(deg2rad*theta_s))
+
+!
+! Direct flux on a flat surface
+!
+  fdir_0 = cste_mars*csza*exp(-taudust/csza)
+
+!
+! Direct flux on the slope
+!
+  fdir = fdir_0 * mu_s/csza
+
+!
+! Reflected flux on the slope
+!
+  fref = albedo * (1-sigma_s) * ftot_0
+
+!
+! Scattered flux on a flat surface
+!
+  fscat_0 = ftot_0 - fdir_0
+
+!
+! Scattering vector (slope vs sky)
+!
+  s_vector=(/ 1., exp(-taudust) , sin(deg2rad*theta_s), sin(deg2rad*theta_s)*exp(-taudust) /)
+
+!
+! Geometry vector (slope vs sun)
+!
+  g_vector=(/ mu_s/csza, 1. /)
+
+!
+! Coupling matrix
+!
+  if (csza .ge. 0.5) then
+	mat_M(:,1) = (/ -0.264,  1.309,  0.208, -0.828 /)
+	mat_M(:,2) = (/  1.291*sigma_s, -1.371*sigma_s, -0.581,  1.641 /)
+	mat_N(:,1) = (/  0.911, -0.777, -0.223,  0.623 /)
+	mat_N(:,2) = (/ -0.933*sigma_s,  0.822*sigma_s,  0.514, -1.195 /)
+
+  else
+	mat_M(:,1) = (/ -0.373,  0.792, -0.095,  0.398 /)
+	mat_M(:,2) = (/  1.389*sigma_s, -0.794*sigma_s, -0.325,  0.183 /)
+	mat_N(:,1) = (/  1.079,  0.275,  0.419, -1.855 /)
+	mat_N(:,2) = (/ -1.076*sigma_s, -0.357*sigma_s, -0.075,  1.844 /)
+  endif
+  !
+  mat_T = mat_M + csza*mat_N
+
+
+!
+! Scattered flux slope ratio
+!
+  if (deg2rad*theta_s <= 0.0872664626) then
+  !
+  ! low angles
+  !
+	s_vector = (/ 1., exp(-taudust) , sin(0.0872664626), sin(0.0872664626)*exp(-taudust) /)
+	ratio = DOT_PRODUCT ( MATMUL( s_vector, mat_T), g_vector )
+        ratio = 1. + (ratio - 1.)*deg2rad*theta_s/0.0872664626
+  else
+  !
+  ! general case
+  !
+	ratio= DOT_PRODUCT ( MATMUL( s_vector, mat_T), g_vector )  
+                  !
+                  ! NB: ratio= DOT_PRODUCT ( s_vector, MATMUL( mat_T, g_vector ) ) is equivalent
+  endif
+
+!
+! Scattered flux on the slope
+!
+  fscat = ratio * fscat_0
+
+endif !! if (csza < 0.01)
+
+!
+! Total flux on the slope
+!
+  ftot = fdir + fref + fscat
+
+!!
+!! Display results
+!!
+!  print *, 'scattered component ', fscat
+!  print *, 'direct component ', fdir
+!  print *, 'reflected component ', fref
+
+END SUBROUTINE param_slope_full
Index: /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/pour_patcher.txt
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/pour_patcher.txt	(revision 51)
+++ /trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/tools_compgcm/pour_patcher.txt	(revision 51)
@@ -0,0 +1,73 @@
+
+****************************************************
+****************************************************
+****************************************************
+****************************************************
+
+
+1. easy part
+diff3 --easy-only --merge oldmeso/meso_physiq.F oldgcm/physiq.F newgcm/physiq.F > meso_physiq.F
+
+2. check quand mm
+vimdiff meso_physiq.F oldmeso/meso_physiq.F
+
+3. reperer conflits
+diff3 --overlap-only meso_physiq.F oldgcm/physiq.F newgcm/physiq.F | more
+--> les lignes sont celles de meso_physiq.F [nouveau fichier]
+[ou] diff3 --show-overlap meso_physiq.F oldgcm/physiq.F newgcm/physiq.F
+
+3bis. reperer conflits
+diff3 --overlap-only --merge meso_physiq.F oldgcm/physiq.F newgcm/physiq.F > dummy
+vimdiff meso_physiq.F dummy
+-- si dummy est OK, alors mv dummy meso_physiq.F
+-- si dummy pas bon, alors n'incorporer que les changements pertinents dans meso_physiq.F
+[ou alors si vous soupconnez l'inverse]
+vimdiff dummy meso_inifis.F
+
+
+4. verifier que les changements sont que mesoscale
+vimdiff meso_physiq.F newgcm/physiq.F
+
+
+
+****************************************************
+****************************************************
+****************************************************
+****************************************************
+****************************************************
+
+
+diff3 --easy-only --merge oldmeso_inifis.F oldgcm_inifis.F newgcm_inifis.F > yeyeye
+diff3 --overlap-only yeyeye oldgcm_inifis.F newgcm_inifis.F
+diff3 --overlap-only --merge yeyeye oldgcm_inifis.F newgcm_inifis.F > yeyeye2
+mv yeyeye2 newmeso_inifis.F
+vimdiff newmeso_inifis.F newgcm_inifis.F [verifier que les changements sont que mesoechelle]
+
+
+
+
+
+
+
+
+
+
+ON VEUT METTRE A JOUR LA PHYSIQUE DU MESO-ECHELLE
+
+1. ancienne version meso
+cp /donnees/aslmd/SVN/trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd_new/libf/phymars/meso_inifis.F .
+
+2. nouvelle version GCM
+cp /u/emlmd/LMDZ.MARS.BETA/libf/phymars/inifis.F .
+
+3. evaluer differences
+\diff -u meso_inifis.F inifis.F > inifis.patch
+
+4. nettoyer inifis.patch en enlevant les changements meso-echelle
+
+5. mettre a jour meso_inifis.F
+patch meso_inifis.F inifis.patch
+
+6. menage 
+\rm inifis.patch inifis.F
+
Index: /trunk/mesoscale/LMD_MM_MARS/prepare_ini
===================================================================
--- /trunk/mesoscale/LMD_MM_MARS/prepare_ini	(revision 50)
+++ /trunk/mesoscale/LMD_MM_MARS/prepare_ini	(revision 51)
@@ -6,4 +6,7 @@
 # Author: Aymeric Spiga - November 2008 #
 #---------------------------------------#
+
+
+echo "If you use new physics, few modifications will be needed"
 
 
Index: /trunk/mesoscale/NOTES.txt
===================================================================
--- /trunk/mesoscale/NOTES.txt	(revision 50)
+++ /trunk/mesoscale/NOTES.txt	(revision 51)
@@ -3,10 +3,11 @@
 - passer aux nouveaux makegcm [en commun avec Ehouarn si on veut le nouveau
   readtesassim qui est en F90]
-- lier gr_fi_dyn qui est dans dyn3d
-- regler le pb du nouveau readtesassim (ou alors le lier tout simplement ou
-  l'appeler meso_readtesassim)
-- regler le pb meso_dustlift (le lier dans makemeso comme point precedent)
 - il faut tester le nest pour verifier les lignes trop longues
 
+(ok) lier gr_fi_dyn qui est dans dyn3d
+(ok) regler le pb du nouveau readtesassim (ou alors le lier tout simplement ou
+  l'appeler meso_readtesassim)
+(ok) regler le pb meso_dustlift (le lier dans makemeso comme point precedent)
+     (car le souci c que dustlift est appele dans vdifc)
 
 RESTE a ADAPTER le LES a la NOUVELLE PHYSIQUE