sxpost.F90

Go to the documentation of this file.

!SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
!SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
!SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
!SFX_LIC for details. version 1.
!     ######spl
        PROGRAM sxpost
!
!!    MODIFICATIONS
!!    -------------
!!      B. Decharme : partition pgd/prep (grid attributes are only in the PGD file)
!!
!-------------------------------------------------------------------------------
!    
USE modd_off_surfex_n
!
        USE modd_io_surf_asc
        USE modd_surf_par
        USE modi_read_surf
        USE modi_get_luout
        USE mode_pos_surf
        USE modd_io_surf_ol, ONLY : xstart,xcount,xstride,lpartr
        
!
        USE yomhook   ,ONLY : lhook,   dr_hook
        USE parkind1  ,ONLY : jprb
!
        USE modi_abor1_sfx
!
        USE modi_get_luout
!
        USE modi_end_io_surf_n
        USE modi_init_io_surf_n
!
        IMPLICIT NONE        

        REAL, ALLOCATABLE, DIMENSION(:)   ::   zloc
        REAL, ALLOCATABLE, DIMENSION(:)   ::   zwrk
        REAL, ALLOCATABLE, DIMENSION(:)   ::   xlon
        REAL, ALLOCATABLE, DIMENSION(:)   ::   xlat
        INTEGER, ALLOCATABLE, DIMENSION(:)::   iwrk2
        CHARACTER(LEN=50)                 ::   ycomment
        CHARACTER(LEN=50)                 ::   nom_article
        CHARACTER(LEN=12)                 ::   hrec
        CHARACTER(LEN=1)                  ::   patchflag
        CHARACTER(LEN=2)                  ::   ypas,ylvl
        CHARACTER(LEN=10)                 ::   cgrid_type
        CHARACTER(LEN=6)                  ::   cmask_save
        LOGICAL                           ::   gfound
        LOGICAL                           ::   linits ! true if PGD has been run
        LOGICAL                           ::   linitp ! true if PREP has been run
        LOGICAL                           ::   lsxnam ! true if SXPOST.nam present
        LOGICAL                           ::   lcoord ! true if LONLAT.dat present
        LOGICAL                           ::   lgeo=.true.  !

        CHARACTER(LEN=28)  :: yluout ='LISTING_SXPOST                '  ! name of listing

        INTEGER    ::   iret
        INTEGER    ::   ini, ini_n
        INTEGER    ::   inj
        INTEGER    ::   if, ic, ip
        INTEGER    ::   ifield, iwfield
        INTEGER    ::   ipatch, jpatch
        INTEGER    ::   ibeg, iend
        REAL(KIND=JPRB) :: zhook_handle

        
        !=====================================================================
        !*
        !** check if file exists
        !*
        !=====================================================================
        IF (lhook) CALL dr_hook('SXPOST',0,zhook_handle)
        CALL surfex_alloc_list(1)

        CALL get_luout('ASCII ',nluout)
        OPEN(unit=nluout,file=adjustl(adjustr(yluout)//'.txt'),&
                        form='FORMATTED',action='WRITE')

        INQUIRE(file='SXPOST.nam', exist=lsxnam)
        IF (.NOT.lsxnam) THEN
           WRITE(*,*)' > SXPOST.nam missing'
           CALL abor1_sfx('SXPOST: NAMELIST SXPOST.nam MISSING')
        ENDIF

        INQUIRE(file='PGD.txt',    exist=linits)
        INQUIRE(file='PREP.txt',   exist=linitp)

        IF (.NOT. linitp .AND. .NOT. linits) THEN
           WRITE(*,*)' NO INPUT FILE FOUND FOR SXPOST'
           WRITE(*,*)' YOU SHOULD AT LEAST RUN PGD! '
           CALL abor1_sfx('SXPOST: NO INPUT FILE')
        ENDIF

        cfilein = 'PGD.txt'

        !=====================================================================
        !*
        !** get number of patches
        !*
        !=====================================================================

         CALL goto_model(1)

 CALL init_io_surf_n(ysurf_cur%DTCO, ysurf_cur%DGU, ysurf_cur%U, &
                        'ASCII ','FULL  ','SURF  ','READ ')          

        CALL read_surf(&
                       'ASCII ','DIM_FULL', ini, iret)
        CALL read_surf(&
                       'ASCII ','GRID_TYPE', cgrid_type, iret)
        CALL read_surf(&
                       'ASCII ','DIM_NATURE', ini_n, iret)        
           
        CALL end_io_surf_n('ASCII ')

        nfull = ini

        IF (ini_n.NE.0) THEN
 CALL init_io_surf_n(ysurf_cur%DTCO, ysurf_cur%DGU, ysurf_cur%U, &
                        'ASCII ','NATURE','SURF  ','READ ')                   

          CALL read_surf(&
                       'ASCII ','PATCH_NUMBER', ipatch, iret)

          CALL end_io_surf_n('ASCII ')
        ENDIF

        !=====================================================================
        !*
        !** get domain size and read latitudes and longitudes 
        !*
        !=====================================================================
 CALL init_io_surf_n(ysurf_cur%DTCO, ysurf_cur%DGU, ysurf_cur%U, &
                        'ASCII ','FULL  ','SURF  ','READ ')          

        OPEN(unit=45,file='SXPOST.nam',form='FORMATTED')
        READ(45,*)ifield

        INQUIRE(file='LONLAT.dat',exist=lcoord)
        ALLOCATE(xlon(ini))
        ALLOCATE(xlat(ini))
        OPEN(unit=30,file='LONLAT.dat',form='FORMATTED')

        IF (lcoord) THEN
          DO ip=1,ini
            READ(30,*)xlon(ip),xlat(ip)
          ENDDO
        ELSE
          WRITE(*,*) 'LONLAT.DAT file missing !'
          IF (cgrid_type=='GAUSS     ') THEN
            CALL posnam(nunit,'FULL  '//' '//'LONGAUSS',gfound,nluout)
          ELSE
            CALL posnam(nunit,'FULL  '//' '//'XLON',gfound,nluout)
          ENDIF

         IF (.NOT.gfound) THEN
           CALL err_stop('XLON            ',cfilein,nluout)
         ELSE
           READ(nunit,fmt=*)
           READ(nunit,fmt='(A50)') ycomment
           READ(nunit,fmt=*,err=100) xlon(:)

           IF (cgrid_type=='GAUSS     ') THEN
             CALL posnam(nunit,'FULL  '//' '//'LATGAUSS',gfound,nluout)
           ELSE
             CALL posnam(nunit,'FULL  '//' '//'XLAT',gfound,nluout)
           ENDIF

           READ(nunit,fmt=*)
           READ(nunit,fmt='(A50)') ycomment
           READ(nunit,fmt=*,err=100) xlat(:)

           DO ip=1,ini
             WRITE(30,*)xlon(ip),xlat(ip)
           ENDDO
         ENDIF
       ENDIF

       CALL end_io_surf_n('ASCII ')

       IF (ifield==0) stop

       !=====================================================================
       !*
       !** read 2d fields from PGD.txt or PREP.txt if exists
       !*
       !=====================================================================

       ! Search var first in PREP file
       ic=0

       DO if=1,ifield

         READ(45,'(A1,1X,A6,1X,A16)') patchflag,cmask,hrec
         cmask_save = cmask

         IF (patchflag == '+') THEN
           inj = ini * ipatch
         ELSE IF (patchflag == '-') THEN
           inj = ini
         ELSE
           print*,' '
           print*,' WRONG PATCHFLAG IN SXPOST.nam '
           print*,' USE + FOR PATCHED VARIABLES '
           print*,' USE - FOR UNPATCHED VARIABLES '
           print*,' '
           print*,' SYNTAX OF SXPOST.nam SHOULD LOOK: '
           print*,' '
           print*,'2 '
           print*,'- FULL   ZS  '
           print*,'+ NATURE TG1 '
           print*,' '
           CALL abor1_sfx('SXPOST: WRONG PATCHFLAG')
         ENDIF

         ALLOCATE(zwrk(inj))
         ic=ic+1
         
         IF (linitp) cfilein = 'PREP.txt'
 CALL init_io_surf_n(ysurf_cur%DTCO, ysurf_cur%DGU, ysurf_cur%U, &
                        'ASCII ',cmask_save,'SURF  ','READ ')           
         CALL posnam(nunit,cmask//' '//hrec,gfound,nluout)
         IF (.NOT.gfound .AND. linitp)THEN
           ! Search now in PGD file
           CALL end_io_surf_n('ASCII ')
           cfilein = 'PGD.txt'
 CALL init_io_surf_n(ysurf_cur%DTCO, ysurf_cur%DGU, ysurf_cur%U, &
                        'ASCII ',cmask_save,'SURF  ','READ ')            
           CALL posnam(nunit,cmask//' '//hrec,gfound,nluout)
         ENDIF
         IF (.NOT.gfound) CALL err_stop(hrec,cfilein,nluout)

         READ(nunit,fmt='(A50)') nom_article
         READ(nunit,fmt='(A50)') ycomment
         READ(nunit,fmt=*,err=100) zwrk
         ALLOCATE(zloc(inj))
         zloc(:)=zwrk(:)
         WHERE(zloc(:)==999.) zloc(:)=-999.

         print*,cmask,' ',hrec,' ','MINVAL = ',minval(zwrk(:)),&
                     ' MAXVAL = ',maxval(zloc(:))
         DEALLOCATE(zloc)

         iwfield=1
         IF (patchflag == '+') THEN
           DO jpatch=1,ipatch
             WRITE(ypas,'(I2)') jpatch
             ylvl=adjustl(ypas(:len_trim(ypas)))
             OPEN(unit=30,file=trim(hrec)//'_p'//trim(ylvl)//'.dat',form='FORMATTED')
             ibeg=ini*(jpatch-1)+1
             iend=ini*jpatch
             IF (lgeo) THEN
               DO ip=1,ini
                 !IF (ZWRK(IP)/=XUNDEF) THEN
                 WRITE(30,*)xlon(ip),xlat(ip),zwrk(ini*(jpatch-1)+ip)
                 !ENDIF
               ENDDO
             ELSE
               WRITE(30,*)nom_article
               WRITE(30,*)ycomment
               DO ip=1,ini
                 WRITE(30,*)zwrk(ip)
               ENDDO
               !WRITE(30,*)INI
               !WRITE(30,'(60F16.8)')ZWRK(IBEG:IEND)
             ENDIF
             CLOSE(30)
           ENDDO
         ELSE
           OPEN(unit=30,file=trim(hrec)//'.dat',form='FORMATTED')
           IF (lgeo) THEN
             DO ip=1,ini
               !IF (ZWRK(IP)/=XUNDEF) THEN
               WRITE(30,*)xlat(ip),xlon(ip),zwrk(ip)
               !ENDIF
             ENDDO
           ELSE
             WRITE(30,*)nom_article
             WRITE(30,*)ycomment
             DO ip=1,ini
               WRITE(30,*)zwrk(ip)
             ENDDO
             !WRITE(30,*)INI
             !WRITE(30,'(60F16.8)')ZWRK(:)
           ENDIF
           CLOSE(30)
         ENDIF

         DEALLOCATE(zwrk)
         CALL end_io_surf_n('ASCII ')

       ENDDO

       CLOSE(nluout)
       CALL surfex_deallo_list
       IF (lhook) CALL dr_hook('SXPOST',1,zhook_handle)


       stop
 100   CONTINUE
       WRITE(nluout,*) ' '
       WRITE(nluout,*) ' ERROR WHEN READING ARTICLE',hrec
       WRITE(nluout,*) ' '
       CLOSE(nluout)
       IF (lhook) CALL dr_hook('SXPOST',1,zhook_handle)

       CONTAINS

       SUBROUTINE err_stop(HREC,CFILEIN,NLUOUT)
       CHARACTER(LEN=12)   ::   hrec
       CHARACTER(LEN=*)    ::   cfilein
       INTEGER             ::   nluout
       REAL(KIND=JPRB) :: zhook_handle
       IF (lhook) CALL dr_hook('ERR_STOP',0,zhook_handle)
       WRITE(nluout,*) ' '
       WRITE(nluout,*) ' ARTICLE ',trim(hrec),' NOT FOUND IN FILE ', cfilein
       WRITE(nluout,*) ' '
       WRITE(*,*) ' '
       WRITE(*,*) ' ARTICLE ',trim(hrec),' NOT FOUND IN FILE ', cfilein
       WRITE(*,*) ' '
       CALL abor1_sfx('SXPOST: ARTICLE '//hrec//' NOT FOUND')
       IF (lhook) CALL dr_hook('ERR_STOP',1,zhook_handle)
       END SUBROUTINE err_stop

       END PROGRAM sxpost