Lire NetCDF.F90

De POLR
Révision datée du 27 juin 2017 à 14:09 par Sennevil (discussion | contributions) (Page créée avec « !*****ISMER-UQAR******************************************************** ! Code du modele ROM en developpement a ISMER-UQAR. ! Lasso ! ! Revision: $Id$ !*******... »)
(diff) ← Version précédente | Voir la version actuelle (diff) | Version suivante → (diff)
Sauter à la navigation Sauter à la recherche

!*****ISMER-UQAR******************************************************** ! Code du modele ROM en developpement a ISMER-UQAR. ! Lasso ! ! Revision: $Id$ !***********************************************************************

     subroutine lire_netcdf(ncid,varid,nt,nj,ni,nk,ki,kf,datanc)

! ! Modification: Simon Senneville - ISMER - 5 octobre 2010 ! Auteur: Simon Senneville - ISMER ! ! Objet: Lecture des fichiers NetCDF ! ! Parametres: ! intent(in): ncid ! Identificateur du fichier NetCDF ! intent(in): varid ! Identificateur de la variable voulue ! intent(in): nt ! dimension relative au temps ! intent(in): nj ! dimension relative a l'axe verticale ! intent(in): ni ! dimension relative a l'axe i ! intent(in): nk ! dimension relative a l'axe k ! intent(in): ki ! position initiale de sur l'axe KSTEP ! intent(in): kf ! position finale de sur l'axe KSTEP ! intent(out): datanc ! champs lu ! ! Retour: n/a ! [ex: 0 : succes] ! [ -1 : erreur] ! ! [note: ] ! :: ! !***********************************************************************

     implicit none

! Data ditionary: declare variables type, definition and units ! type :: var_name ! definition, units ! [type,dimension([x,y]) :: var_name ! definition, units] ! [type,intent(in|out|inout]):: name ! definition, units] ...

     integer,intent(in)     :: nt, ni, nj, nk, ki, kf, varid, ncid
     real,dimension(nt,nj,ni,nk),intent(out) :: datanc
     real,dimension(nk,ni)  :: datanc_no_comp
     integer                :: i, k

!-----------------------------------------------------------------------

     integer                :: idims, it, iw
     integer                :: ret, xtype, ndims, natts
     integer                :: id, wi2d, wi3d, dimw, iwid, indx, indy, indz, indice
     integer                :: dimx, dimy, dimz, dimt, stat, dimxy
     integer                :: xid, yid, zid, tid, tvid, ksvid, dtid, dimtid, dimwid
     integer                :: tread_ti, tread_tf
     integer                :: deltat_nc, dt_in,kstart, kstart_p1, dkstep
     integer,dimension(2)   :: indstart_cp, indcount_cp
     integer,dimension(6)   :: ncref                      ! Date de reference (ss,mm,hh,j,m,a)
     real(kind=8)           :: nct, ncstart, ncend, ncvalid
     real,allocatable,dimension(:,:)  :: data_wet
     integer,allocatable,dimension(:) :: dimids
     integer,allocatable,dimension(:) :: wet
     character(len=80)      :: nom 
     character(len=80),allocatable,dimension(:) :: dimnames
     logical                :: v2d, v3d, comp, vtd
                                                          ! Including the NetCDF library for fortran
  1. include "netcdf.inc"
     vtd  = .FALSE.                                       ! Initialisation de la semaphore vtd
     ret = nf_inq_varndims(ncid,varid,ndims)              ! Trouver le nombre de dimension de la variable TAG
     if (ret.ne.NF_NOERR) then
        write(*,*) nf_strerror(ret)
        stop
     end if
     allocate (dimids(ndims),dimnames(ndims),stat=stat)   ! Allouer la dimension de dimids et dimnames
     if (stat /= 0) then
        write(*,*)'Erreur: allocate'
        stop
     end if
     ret = nf_inq_var(ncid,varid,nom,xtype,ndims,dimids(1),natts) ! Lire les informaitons relatives a la variable varid
     if (ret.ne.NF_NOERR) then
        write(*,*) nf_strerror(ret)
        stop
     end if
     do idims = 1,ndims
       ret = nf_inq_dimname(ncid,dimids(idims),nom)       ! Lire de nom des dimensions
       if (ret.ne.NF_NOERR) then
          write(*,*) nf_strerror(ret)
          stop
       end if
       dimnames(idims) = nom(:)
       if (trim(dimnames(idims))=='NVALID2D') then        ! Verifier si c'est une variable 2D compressee
         comp = .TRUE.
         v2d  = .TRUE.
         ret = nf_inq_dimid(ncid,'NVALID2D',dimwid)       ! Si oui, lire l'id de NVALID2D
         if (ret.ne.NF_NOERR) then
            write(*,*) nf_strerror(ret)
            stop
         end if
         ret = nf_inq_dimlen(ncid,dimwid,dimw)            ! lire la dimension de NVALID2D
         if (ret.ne.NF_NOERR) then
            write(*,*) nf_strerror(ret)
            stop
         end if
         ret = nf_inq_varid(ncid,'INDEXVALID2D',iwid)     ! lire l'id de INDEXVALID2D
         if (ret.ne.NF_NOERR) then
            write(*,*) nf_strerror(ret)
            stop
         end if
       end if
       if (trim(dimnames(idims))=='NVALID3D') then        ! Verifier si c'est une variable 3D compressee
         comp = .TRUE.
         v3d  = .TRUE.
         ret = nf_inq_dimid(ncid,'NVALID3D',dimwid)       ! Si oui, lire l'id de NVALID3D
         if (ret.ne.NF_NOERR) then
            write(*,*) nf_strerror(ret)
            stop
         end if
         ret = nf_inq_dimlen(ncid,dimwid,dimw)            ! lire la dimension de NVALID3D
         if (ret.ne.NF_NOERR) then
            write(*,*) nf_strerror(ret)
            stop
         end if
         ret = nf_inq_varid(ncid,'INDEXVALID3D',iwid)     ! lire l'id de INDEXVALID3D
         if (ret.ne.NF_NOERR) then
            write(*,*) nf_strerror(ret)
            stop
         end if
       end if
       if (trim(dimnames(idims))=='TROM') then            ! Verifier si c'est une variable fonction de TROM
         vtd  = .TRUE.
         ret = nf_inq_dimid(ncid,'TROM',tid)              ! Si oui, lire l'id de TROM
         if (ret.ne.NF_NOERR) then
           write(*,*) nf_strerror(ret)
           stop
         end if
         ret = nf_inq_varid(ncid,'TROM',dimt)             ! lire l'id variable TROM 
         if (ret.ne.NF_NOERR) then
           write(*,*) nf_strerror(ret)
           stop
         end if
       end if
       if (trim(dimnames(idims))=='KROM') then            ! Non-compressee
         comp  = .FALSE.
       end if
       if (trim(dimnames(idims))=='IROM') then            ! Non-compressee
         comp  = .FALSE.
       end if
       if (trim(dimnames(idims))=='ZROM') then            ! Non-compressee
         comp  = .FALSE.
       end if
     end do
     if (comp) then                                       ! Si la variable est compressee, il faut la decompresse
       allocate (wet(dimw),stat=stat)                     ! Alloue l'espace a l'axe de comrpession
       if (stat/=0) then
         write(*,*) 'Error with allocate of the variable: wet'
         stop
       end if
       allocate (data_wet(dimw,nt),stat=stat)             ! Alloue l'espace a la variable compressee
       if (stat/=0) then
         write(*,*) 'Error with allocate of the variable: data_wet'
         stop
       end if
       ret = nf_get_var_int(ncid,iwid,wet)                ! Lire la variable de l'axe de compression
       if (ret.ne.NF_NOERR) then
          write(*,*) nf_strerror(ret)
          stop
       end if
                                                          !------------------------------
       if (vtd) then                                      ! Si la variable est en fonction de TROM
         indstart_cp(1) = 1                               ! Lire une ecriture
         indstart_cp(2) = ki                              ! Le ki ieme element
         indcount_cp(1) = dimw                            ! de dimw elements 
         indcount_cp(2) = nt                              ! pour nt pas de temps
         ret = nf_get_vara_real(ncid,varid,indstart_cp(1),indcount_cp(1), data_wet(1,1)); ! Lire cette information dans data_wet
         if (ret.ne.NF_NOERR) then
            write(*,*) nf_strerror(ret)
            stop
         end if
                                                          !------------------------------
       else                                               ! Si la variable est independante du temps
         indstart_cp(1) = 1                               ! Lire une ecriture
         indstart_cp(2) = 1                               ! Le ki ieme element
         indcount_cp(1) = dimw                            ! de dimw elements 
         indcount_cp(2) = 1                               ! pour nt pas de temps
         ret = nf_get_vara_real(ncid,varid,indstart_cp(1),indcount_cp(1), data_wet(1,1)); ! Lire cette information dans data_wet
         if (ret.ne.NF_NOERR) then
            write(*,*) nf_strerror(ret)
            stop
         end if
       end if
                                                          !------------------------------
       dimxy=ni*nk;                                       ! Dimension horizontale
       do it = 1 , nt                                     ! Decompression du champs compressee
         do iw = 1 , dimw
           indice = wet(iw);
           indz   = (indice-1) / dimxy;
           indice = indice-indz*dimxy;
           indy   = (indice-1) / nk;
           indx   = indice-indy*nk;
           indy   = indy + 1 ;
           indz   = indz + 1 ;
           datanc(it,indz,indy,indx) = data_wet(iw,it);   ! La variable decompresse est lu dans datanc
         end do
       end do
       deallocate(wet,data_wet)
     else                                                 ! Si la variable ne depend pas du temps
       ret = nf_get_var_real(ncid,varid,datanc_no_comp)   ! Lire cette information dans datanc directment
                                                          ! nf_get_var_real lit les variable dans un ordre different (comme en C)
       if (ret.ne.NF_NOERR) then
          write(*,*) nf_strerror(ret)
          stop
       endif
       forall (i=1:ni,k=1:nk)                             ! Transformation de (k,i) a (i,k)
         datanc(1,1,i,k)=datanc_no_comp(k,i)              ! A_FAIRE recoder cette partie afin qu'elle soit plus robuste
       endforall
     endif
     ret = nf_close(ncid)                                 ! Fermeture du fichier ncid ouvert par getdims
     if (ret.ne.NF_NOERR) then
        write(*,*) nf_strerror(ret)
        stop
     end if
     deallocate(dimids,dimnames,stat=stat)                ! Deallocation de l'espace memoire
     if (stat /= 0) then
        write(*,*)'Dans lire_netcdf Erreur: deallocate'
        stop
     end if
                                                          !------------------------------
     end subroutine lire_netcdf