source: trunk/LMDZ.MARS/libf/phymars/NLTEdlvr09_ZGRID.F @ 414

c***********************************************************************

        subroutine NLTEdlvr09_ZGRID (n_gcm,  
     @          p_gcm, t_gcm, z_gcm, 
     @          co2vmr_gcm, n2vmr_gcm, covmr_gcm, o3pvmr_gcm ,mmean_gcm,
     @          cpnew_gcm)

c       jul 2011 malv+fgg    First version
c***********************************************************************
                                                
        implicit none                                  
        
        include "dimensions.h"
        include "dimphys.h"
        include 'nltedefs.h'          
        include 'nlte_atm.h'        
        include 'tcr_15um.h'
        include 'nlte_data.h'       
        include 'chimiedata.h'
        include 'conc.h'
                                                
c Arguments 
        integer n_gcm
        real p_gcm(n_gcm), t_gcm(n_gcm)
        real co2vmr_gcm(n_gcm), n2vmr_gcm(n_gcm)
        real covmr_gcm(n_gcm), o3pvmr_gcm(n_gcm)
        real z_gcm(n_gcm)
        real mmean_gcm(n_gcm)
        real cpnew_gcm(n_gcm)

c local variables                               
        integer i, j  , iz                               
!       real  distancia, meanm, gz, Hkm
        real  zmin, zmax, deltazz, deltazzy
        real  nt_gcm(n_gcm)
        real  mmean_nlte(n_gcm),cpnew_nlte(n_gcm)
              
c functions                                     
        external        hrkday_convert                       
        real            hrkday_convert                          
                                                
c***********************************************************************


! Define working grid for MZ1D model (NL, ZL, ZMIN) 
! y otro mas fino para M.Curtis (NZ, ZX, ZXMIN = ZMIN

! Para ello hace falta una z de ref del GCM, que voy a suponer la inferior

! Primero, construimos escala z_gcm 

!       z_gcm (1) = zmin_gcm             ! [km]

        !write (*,*) ' iz, p, g, H, z =', 1, p_gcm(1), z_gcm(1)
!       do iz = 2, n_gcm
!       do iz=1,n_gcm
!          z_gcm(iz)=zlay(iz)/1.e3

!         meanm = ( co2vmr_gcm(iz)*44. + o3pvmr_gcm(iz)*16. 
!     @               + n2vmr_gcm(iz)*28. + covmr_gcm(iz)*28. ) 
!         meanm = meanm / n_avog
!         distancia = ( radio + z_gcm(iz-1) )*1.e5 
!         gz = gg * masa / ( distancia * distancia ) 
!          Hkm = 0.5*( t_gcm(iz)+t_gcm(iz-1) ) / ( meanm * gz )
!          Hkm = kboltzman * Hkm *1e-5                           ! [km] 
!          z_gcm(iz) = z_gcm(iz-1) - Hkm * log( p_gcm(iz)/p_gcm(iz-1) )

          !write (*,*) iz, p_gcm(iz), gz, Hkm, z_gcm(iz)

!        enddo
! Segundo, definimos los límites del modelo, entre las 2 presiones clave

        ! Bottom boundary for NLTE model : Pbottom=2e-2mb=1.974e-5 atm
        jlowerboundary = 1 
        do while ( p_gcm(jlowerboundary) .gt. Pbottom_atm )
           jlowerboundary = jlowerboundary + 1
        enddo
        zmin = z_gcm(jlowerboundary)
!       write (*,*) ' jlowerboundary, Pmin, zmin =', 
!     @            jlowerboundary, p_gcm(jlowerboundary), zmin

        ! Top boundary for NLTE model : Ptop=2e-7mb = 1.974e-5 atm
        jtopboundary = jlowerboundary  
        do while ( p_gcm(jtopboundary) .gt. Ptop_atm ) 
           jtopboundary = jtopboundary + 1
        enddo
        zmax = z_gcm(jtopboundary)
!       write (*,*) ' jtopboundary, Pmax, zmax =', 
!     @      jtopboundary, p_gcm(jtopboundary),zmax

        deltaz = (zmax-zmin) / (nl-1) 
        do i=1,nl                                      
            zl(i) = zmin + (i-1) * deltaz
        enddo                                          
!       write (*,*) ' ZL grid:  dz,zmin,zmax ', deltaz, zl(1),zl(nl)
! Creamos el perfil interpolando
        call intersp (    pl,zl,nl,      p_gcm,z_gcm,n_gcm, 2)   ! [atm]
        call intersp (     t,zl,nl,      t_gcm,z_gcm,n_gcm, 1)       
        do i = 1, n_gcm
           nt_gcm(i) = 7.339e+21 * p_gcm(i) / t_gcm(i)    ! [cm-3]         
        enddo
        call intersp (    nt,zl,nl,     nt_gcm,z_gcm,n_gcm, 2)       
        call intersp (co2vmr,zl,nl, co2vmr_gcm,z_gcm,n_gcm, 1)       
        call intersp ( n2vmr,zl,nl,  n2vmr_gcm,z_gcm,n_gcm, 1)       
        call intersp ( covmr,zl,nl,  covmr_gcm,z_gcm,n_gcm, 1)       
        call intersp (o3pvmr,zl,nl, o3pvmr_gcm,z_gcm,n_gcm, 1)  
        call intersp (mmean_nlte,zl,nl,mmean_gcm,z_gcm,n_gcm,1)
        call intersp (cpnew_nlte,zl,nl,cpnew_gcm,z_gcm,n_gcm,1)
        

        do i = 1, nl

           co2(i) = nt(i) * co2vmr(i)
           n2(i) = nt(i) * n2vmr(i)
           co(i) = nt(i) * covmr(i)
           o3p(i) = nt(i) * o3pvmr(i)

!               hrkday_factor(i) =  hrkday_convert( t(i),        
!     @           co2vmr(i), o3pvmr(i), n2vmr(i), covmr(i) )
           hrkday_factor(i) = hrkday_convert (mmean_nlte(i),cpnew_nlte(i))

        enddo
                                                
                                              

c  Fine grid for transmittance calculations

        deltazy = (zmax-zmin) / (nzy-1) 
        do i=1,nzy                                      
            zy(i) = zmin + (i-1) * deltazy                   
        enddo                                          
!       write (*,*) ' ZY grid:  nzy,dzy,zmin,zmax ', 
!     @         nzy, deltazy, zy(1),zy(nzy)

        call intersp (    py,zy,nzy,      p_gcm,z_gcm,n_gcm, 2)   ! [atm]
        call intersp (    ty,zy,nzy,      t_gcm,z_gcm,n_gcm, 1)       
        call intersp (   nty,zy,nzy,     nt_gcm,z_gcm,n_gcm, 2)       

        call intersp (  co2y,zy,nzy,   co2vmr_gcm,z_gcm,n_gcm, 1) 
        do i=1,nzy 
           co2y(i) = co2y(i) * nty(i)
        enddo


  

c end                                           
        return                                         
        end