Subversion Repositories pvinversion.ecmwf

      PROGRAM Cutnetcdf

c     -----------------------------------------------------------------------
c     Cut and split an input netcdf 
c     Michael Sprenger / Summer 2006
c     -----------------------------------------------------------------------

      implicit none

c     -----------------------------------------------------------------------
c     Declaration of variables
c     -----------------------------------------------------------------------

c     Maximum domain size
      integer        nxmax,nymax,nzmax,ntmax
      parameter      (nxmax=400,nymax=400,nzmax=400,ntmax=5)

c     Variables for input and output netcdf file
      integer        cdfid
      real           phymin(4),phymax(4),stag(4), misdat
      integer        ndim, vardim(4),  error, nlev
      integer        cstid
      real           dx, dy
      real           aklev(nzmax),bklev(nzmax),aklay(nzmax),bklay(nzmax)
      real           pollon, pollat
      integer        stdate(5)
      integer        ncdfid
      integer        datar(14)
      character*80   cfn
      real           varold(nxmax*nymax*nzmax)
      real           time(ntmax)

c     Parameters
      character*80   oldfile(400),newfile(400)
      character*80   oldname(400),newname(400)
      integer        novars

c     Auxiliary variables
      integer        i,j,k,l,t,ntimes
      integer        imin, imax, jmin, jmax, tmin, tmax
      integer        nx,ny,nz
      integer        crename(400),crefile(400)
      integer        isnew
      character*80   varname

c     -----------------------------------------------------------------------
c     Preparations
c     -----------------------------------------------------------------------

      print*,'*********************************************************'
      print*,'* cutnetcdf                                             *'
      print*,'*********************************************************'

c     Read entries from argument file
      open(10,file='fort.10')
        novars=1
 100    read(10,*,end=110) oldname(novars),newname(novars),
     >                     oldfile(novars),newfile(novars)
        novars=novars+1
        goto 100
 110  close(10)
      novars=novars-1

c     Init the flags for creating files and variables
      do i=1,novars
         crename(i)=0
         crefile(i)=0
      enddo

c     -----------------------------------------------------------------------
c     Loop through data points
c     -----------------------------------------------------------------------

      do t=1,novars  

c       Write info 
        write(*,'(a10,a5,3x,a5,3x,a15,3x,a15)')
     >           '  Split  ',trim(oldname(t)),
     >                       trim(newname(t)),
     >                       trim(oldfile(t)),
     >                       trim(newfile(t))

c       Open input file and read some parametersq
        call cdfopn(oldfile(t),cdfid,error)
        if (error.ne.0) goto 997
        call getcfn(cdfid,cfn,error)
        if (error.ne.0) goto 997
        call cdfopn(cfn,cstid,error)
        if (error.ne.0) goto 997
        call gettimes(cdfid, time, ntimes, error) 
        if (error.ne.0) goto 997
        call getgrid(cstid, dx, dy, error)
        if (error.ne.0) goto 997
        call getlevs(cstid, nlev, aklev, bklev, aklay, bklay, error)
        if (error.ne.0) goto 997
        call getpole(cstid,pollon,pollat,error)
        if (error.ne.0) goto 997        
        call getstart(cstid,stdate,error)
        if (error.ne.0) goto 997

c       Set new grid parameters and read data of subdomain
        call getdef(cdfid,oldname(t),ndim,misdat,
     >              vardim,phymin,phymax,stag,error)
        if (error.ne.0) goto 997
        nx=vardim(1)
        ny=vardim(2)
        nz=vardim(3)  

c       Load data
        varname=oldname(t)
        call getdat(cdfid,varname,time(1),0,varold,error)
        if (error.ne.0) goto 997

c       Create file if necessary (otherwise open it for writing)
        isnew=1
        do k=1,novars
           if ((crefile(k).eq.1).and.
     >         (newfile(t).eq.newfile(k))) then
              isnew=0
           endif
        enddo
        if (isnew.eq.1) then
           cfn=trim(newfile(t))//'_cst'
           datar(1)  = vardim(1)
           datar(2)  = vardim(2)
           datar(3)  = 1000.*phymax(2)
           datar(4)  = 1000.*phymin(1)
           datar(5)  = 1000.*phymin(2)
           datar(6)  = 1000.*phymax(1)
           datar(7)  = 1000.*dx
           datar(8)  = 1000.*dy
           datar(9)  = nlev
           datar(10) = 1
           datar(11) = 0
           datar(12) = 0
           datar(13) = 1000.*pollon
           datar(14) = 1000.*pollat
           call wricst(cfn,datar,
     >                 aklev,bklev,aklay,bklay,stdate)
           call crecdf(trim(newfile(t)),ncdfid, phymin, phymax, 
     >                 ndim,cfn,error)
           if (error.ne.0) goto 998
        else
           call cdfwopn(trim(newfile(t)),ncdfid,error)
           if (error.ne.0) goto 998
        endif 
        do k=1,novars
           if (newfile(k).eq.newfile(t)) then
              crefile(k)=1
           endif
        enddo

c       Create variable if necessary
        isnew=1
        do k=1,novars
           if ((newfile(t).eq.newfile(k)).and.
     >         (newname(t).eq.newname(k)).and.
     >         (crename(t).eq.1)) then
              isnew=0
           endif
        enddo
        if (isnew.eq.1) then
           call putdef(ncdfid, newname(t), ndim, misdat, vardim, 
     >                  phymin, phymax, stag, error) 
           if (error.ne.0) goto 998
        endif
        do k=1,novars
           if ((newname(k).eq.newname(t)).and.
     >         (crefile(k).eq.1)) then
              crename(k)=1
           endif
        enddo

c       Write data
        call putdat(ncdfid,newname(t),time(1),0,varold,error)
        if (error.ne.0) goto 998

c       Close files
        call clscdf(cdfid, error)
        if (error.ne.0) goto 997
        call clscdf(ncdfid, error)
        if (error.ne.0) goto 998
        call clscdf(cstid, error)
        if (error.ne.0) goto 997

      enddo

c     -----------------------------------------------------------------------
c     Exception handling
c     -----------------------------------------------------------------------

      stop

 997  print*,'Problems with input file... Stop'
      stop

 998  print*,'Problems with output file... Stop'
      stop


      end