c&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
c Modele approche de rayonnement: methode pertubationnelle
c dst 1er ordre en kgaz et ksuies
c --------------------------------
c Amaury janv 2001
c&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
c
c
c.......................................................................

      program mosim


      implicit none
      include 'cecile.inc'
      include 'propradia.inc'
      include 'entre.inc'
      integer itmp1, i,j,m
      double precision bil ,var_bil, psid, supp
      double precision fvnew, fmnew, tempnew

      dimension fvnew(1:n_mx)
      dimension tempnew(0:n_mx+1)
      dimension fmnew(1:ngaz_mx, 0:n_mx+1)

      dimension bil(0:n_mx+1), var_bil(0:n_mx+1), psid(0:n_mx,0:n_mx)
      dimension supp(0:n_mx+1,0:n_mx+1)

      character*32 chemin

      print *, 'CHANGER LE CHEMIN DE LA SOL DE REFERENCE'
      read *
      chemin = '../raycecile/resultats/mcmh2o_05'
      write(*,*) 'n_mx=',n_mx 
      do i=0,n+1
         do j=0,n+1
               supp(i,j)=0.D+0
         enddo
      enddo
c ************************************************************************
      print *, 'indice'
      open(20,file=chemin // '/ind.out', status
     & ='unknown')
      read(20,*) n, rsup
      read(20,*) n, rsup
      read(20,*) n, rsup
      close(20)
c...............................................................................
      print *, 'psi'
      open(20,file=chemin // '/psi.out',status
     & ='unknown')
      do i = 0,n+1
              do j=0,n+1
                 read(20,*) in, iin, r(i), r(j)
     & ,psi(i,j), var_psi(i,j)
              enddo
              read(20,*)
      enddo
      close (20)

c...............................................................................
      print *, 'dpsi'
      open(20,file=chemin // '/dpsi.out',status
     & ='unknown')
      do i=0,n+1
                do j=0,n+1
                        do m=1,n
                                read(20,'(3(i4),4(e22.13))') in, iin,
     &  itmp1, dpsi6dfv(i,j,m)
                        enddo
                        read(20,*)
                enddo
      enddo
      close(20)
c...............................................................................
      print *, 'dpsidfm'
       open(20,file=chemin // '/dpsidfm.out',
     & status='unknown')
        do in=0,n+1
                do iin=0,n+1
                        do itmp1=1,n
                                read(20,'(6(e22.13))')
     &   dpsi6dfmco(in,iin,itmp1),
     &   dpsi6dfmco2(in,iin,itmp1),
     &   dpsi6dfmh2o(in,iin,itmp1)
                        enddo
                enddo
      enddo
       close(20)

c...............................................................................
      print *, 'dpidt'
      open(20,file=chemin // '/dpsidt.out',status
     & ='unknown')
      do i=0,n+1
                do j=0,n+1
                        do m=1,n
                                read(20,*) in, iin,
     &   itmp1, dpsi6dt(i,j,m)
                        enddo
                        read(20,*)
                enddo
      enddo
      close(20)
c-----------------------------------------------------------------
      print *, 'fv'
       open(20,file=chemin // '/fv.in',status='unknown')
       do i=1,n
          read (20,*) fv(i)
       enddo
       close(20)
      print *, 'fmco'
        open(20,file=chemin // '/fmco.in',
     &     status='unknown')
       do i=1,n
          read (20,*) fm(1,i)
       enddo
       close(20)
      print *, 'fmco2'
        open(20,file=chemin // '/fmco2.in',
     &      status='unknown')
       do i=1,n
          read (20,*) fm(2,i)
       enddo
       close(20)
      print *, 'fmh2o'
        open(20,file=chemin // '/fmh2o.in',
     &      status='unknown')
       do i=1,n
          read (20,*) fm(3,i)
       enddo
       close(20)
      print *, 'temp'
       open(20,file=chemin // '/temperat.in',
     &     status='unknown')
       do i=0,n+1
          read (20,*) temp(i)
       enddo

       close(20)
c---------------------------------------------------------------------
      print *, 'NEW IN fv'
       open(20,file='./fv.in',status='unknown')
       do i=1,n
          read (20,*) fvnew(i)
       enddo
       close(20)
      print *, 'NEW IN fmco'
        open(20,file='./fmco.in',status='unknown')
       do i=1,n
          read (20,*) fmnew(1,i)
       enddo
       close(20)
      print *, 'NEW IN fmco2'
        open(20,file='./fmco2.in',status='unknown')
       do i=1,n
          read (20,*) fmnew(2,i)
       enddo
       close(20)
      print *, 'NEW IN fmh2o'
        open(20,file='./fmh2o.in',status='unknown')
       do i=1,n
          read (20,*) fmnew(3,i)
       enddo
       close(20)
      print *, 'NEW IN temp'
       open(20,file='./temperat.in',status='unknown')
       do i=0,n+1
          read (20,*) tempnew(i)
       enddo

       close(20)

       print *, 'FIN LECTURE-goto DL'
       do i=0,n+1
          do j=0,n+1
             do m=1,n
                supp(i,j)=supp(i,j)+dpsi6dfv(i,j,m)*(fvnew(m)- fv(m))
     &                         +dpsi6dfmco(i,j,m)*(fmnew(1,m)-fm(1,m))
     &                         +dpsi6dfmco2(i,j,m)*(fmnew(2,m)-fm(2,m))
     &                         +dpsi6dfmh2o(i,j,m)*(fmnew(3,m)-fm(3,m))
     &                         +dpsi6dt(i,j,m)*(tempnew(m)-temp(m))
             enddo
             psi(i,j)= psi(i,j)+supp(i,j)                      
          enddo
      enddo
      call fich

      end


