Changeset b4d29ce in flexpart.git for src/readoutgrid.f90

Timestamp:

Jul 3, 2014, 2:55:50 PM (10 years ago)

Author:

Harald Sodemann <harald.sodemann@…>

Branches:

master, 10.4.1_pesei, FPv9.3.1, FPv9.3.1b_testing, FPv9.3.2, GFS_025, NetCDF, bugfixes+enhancements, deposition, dev, fp9.3.1-20161214-nc4, grib2nc4_repair, inputlist, laptop, release-10, release-10.4.1, scaling-bug, svn-petra, svn-trunk, univie

Children:

4b093cc, 75dfd42, 326b31b

Parents:

87910af

Message:

• Implemented optional namelist input for COMMAND, RELEASES, SPECIES, AGECLASSES,OUTGRID,OUTGRID_NEST,RECEPTORS
• Implemented com_mod switch nmlout to write input files as namelist to the output directory (.true. by default)
• Proposed updated startup and runtime output (may change back to previous info if desired)

git-svn-id: ​http://flexpart.flexpart.eu:8088/svn/FlexPart90/trunk@27 ef8cc7e1-21b7-489e-abab-c1baa636049d

File:

• src/readoutgrid.f90 (modified) (6 diffs)

src/readoutgrid.f90

@@ -29 +29 @@
   !                                                                            *
   !     4 June 1996                                                            *
+  !     HSO, 1 July 2014
+  !     Added optional namelist input
   !                                                                            *
   !*****************************************************************************
@@ -53 +55 @@
   real,parameter :: eps=1.e-4
 
-
+  ! namelist variables
+  integer, parameter :: maxoutlev=500
+  integer :: readerror
+  real,allocatable, dimension (:) :: outheights
+
+  ! declare namelist
+  namelist /outgrid/ &
+    outlon0,outlat0, &
+    numxgrid,numygrid, &
+    dxout,dyout, &
+    outheights
+
+  ! allocate large array for reading input
+  allocate(outheights(maxoutlev),stat=stat)
+  if (stat.ne.0) write(*,*)'ERROR: could not allocate outheights'
+
+  ! helps identifying failed namelist input
+  dxout=-1.0
+  outheights=-1.0
 
   ! Open the OUTGRID file and read output grid specifications
   !**********************************************************
 
-  open(unitoutgrid,file=path(1)(1:length(1))//'OUTGRID',status='old', &
-       err=999)
-
-
-  call skplin(5,unitoutgrid)
-
-
-  ! 1.  Read horizontal grid specifications
-  !****************************************
-
-  call skplin(3,unitoutgrid)
-  read(unitoutgrid,'(4x,f11.4)') outlon0
-  call skplin(3,unitoutgrid)
-  read(unitoutgrid,'(4x,f11.4)') outlat0
-  call skplin(3,unitoutgrid)
-  read(unitoutgrid,'(4x,i5)') numxgrid
-  call skplin(3,unitoutgrid)
-  read(unitoutgrid,'(4x,i5)') numygrid
-  call skplin(3,unitoutgrid)
-  read(unitoutgrid,'(4x,f12.5)') dxout
-  call skplin(3,unitoutgrid)
-  read(unitoutgrid,'(4x,f12.5)') dyout
-
+  open(unitoutgrid,file=path(1)(1:length(1))//'OUTGRID',status='old',form='formatted',err=999)
+
+  ! try namelist input
+  read(unitoutgrid,outgrid,iostat=readerror)
+  close(unitoutgrid)
+
+  if ((dxout.le.0).or.(readerror.ne.0)) then
+
+    readerror=1
+
+    open(unitoutgrid,file=path(1)(1:length(1))//'OUTGRID',status='old',err=999)
+
+    call skplin(5,unitoutgrid)
+
+    ! 1.  Read horizontal grid specifications
+    !****************************************
+
+    call skplin(3,unitoutgrid)
+    read(unitoutgrid,'(4x,f11.4)') outlon0
+    call skplin(3,unitoutgrid)
+    read(unitoutgrid,'(4x,f11.4)') outlat0
+    call skplin(3,unitoutgrid)
+    read(unitoutgrid,'(4x,i5)') numxgrid
+    call skplin(3,unitoutgrid)
+    read(unitoutgrid,'(4x,i5)') numygrid
+    call skplin(3,unitoutgrid)
+    read(unitoutgrid,'(4x,f12.5)') dxout
+    call skplin(3,unitoutgrid)
+    read(unitoutgrid,'(4x,f12.5)') dyout
+
+  endif
 
   ! Check validity of output grid (shall be within model domain)
@@ -91 +119 @@
   if ((outlon0+eps.lt.xlon0).or.(outlat0+eps.lt.ylat0) &
        .or.(xr.gt.xr1+eps).or.(yr.gt.yr1+eps)) then
-    write(*,*) 'outlon0,outlat0:'
     write(*,*) outlon0,outlat0
-    write(*,*) 'xr1,yr1,xlon0,ylat0,xr,yr,dxout,dyout:'
     write(*,*) xr1,yr1,xlon0,ylat0,xr,yr,dxout,dyout
     write(*,*) ' #### FLEXPART MODEL ERROR! PART OF OUTPUT    ####'
@@ -104 +130 @@
   ! 2. Count Vertical levels of output grid
   !****************************************
-  j=0
-100   j=j+1
+
+  if (readerror.ne.0) then
+    j=0
+100 j=j+1
     do i=1,3
       read(unitoutgrid,*,end=99)
@@ -112 +140 @@
     if (outhelp.eq.0.) goto 99
     goto 100
-99   numzgrid=j-1
-
-    allocate(outheight(numzgrid) &
-         ,stat=stat)
-    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
-    allocate(outheighthalf(numzgrid) &
-         ,stat=stat)
-    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
-
-
-  rewind(unitoutgrid)
-  call skplin(29,unitoutgrid)
+99  numzgrid=j-1
+  else
+    do i=1,maxoutlev
+      if (outheights(i).lt.0) exit
+    end do
+    numzgrid=i-1
+  end if
+
+  allocate(outheight(numzgrid),stat=stat)
+  if (stat.ne.0) write(*,*)'ERROR: could not allocate outheight'
+  allocate(outheighthalf(numzgrid),stat=stat)
+  if (stat.ne.0) write(*,*)'ERROR: could not allocate outheighthalf'
 
   ! 2. Vertical levels of output grid
   !**********************************
 
-  j=0
-1000   j=j+1
-    do i=1,3
-      read(unitoutgrid,*,end=990)
-    end do
-    read(unitoutgrid,'(4x,f7.1)',end=990) outhelp
-    if (outhelp.eq.0.) goto 99
-    outheight(j)=outhelp
-    goto 1000
-990   numzgrid=j-1
-
+  if (readerror.ne.0) then
+
+    rewind(unitoutgrid)
+    call skplin(29,unitoutgrid)
+
+    do j=1,numzgrid
+      do i=1,3
+        read(unitoutgrid,*)
+      end do
+      read(unitoutgrid,'(4x,f7.1)') outhelp
+      outheight(j)=outhelp
+      outheights(j)=outhelp
+    end do
+    close(unitoutgrid)
+
+  else
+
+    do j=1,numzgrid
+      outheight(j)=outheights(j)
+    end do
+
+  endif
+
+  ! write outgrid file in namelist format to output directory if requested
+  if (nmlout.eqv..true.) then
+    ! reallocate outheights with actually required dimension for namelist writing
+    deallocate(outheights)
+    allocate(outheights(numzgrid),stat=stat)
+    if (stat.ne.0) write(*,*)'ERROR: could not allocate outheights'
+
+    do j=1,numzgrid
+      outheights(j)=outheight(j)
+    end do
+
+    open(unitoutgrid,file=path(2)(1:length(2))//'OUTGRID.namelist',err=1000)
+    write(unitoutgrid,nml=outgrid)
+    close(unitoutgrid)
+  endif
 
   ! Check whether vertical levels are specified in ascending order
@@ -160 +215 @@
   end do
 
-
   xoutshift=xlon0-outlon0
   youtshift=ylat0-outlat0
-  close(unitoutgrid)
-
-    allocate(oroout(0:numxgrid-1,0:numygrid-1) &
-         ,stat=stat)
-    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
-    allocate(area(0:numxgrid-1,0:numygrid-1) &
-         ,stat=stat)
-    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
-    allocate(volume(0:numxgrid-1,0:numygrid-1,numzgrid) &
-         ,stat=stat)
-    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
-    allocate(areaeast(0:numxgrid-1,0:numygrid-1,numzgrid) &
-         ,stat=stat)
-    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
-    allocate(areanorth(0:numxgrid-1,0:numygrid-1,numzgrid) &
-         ,stat=stat)
-    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
+
+  allocate(oroout(0:numxgrid-1,0:numygrid-1),stat=stat)
+  if (stat.ne.0) write(*,*)'ERROR: could not allocate oroout'
+  allocate(area(0:numxgrid-1,0:numygrid-1),stat=stat)
+  if (stat.ne.0) write(*,*)'ERROR: could not allocate area'
+  allocate(volume(0:numxgrid-1,0:numygrid-1,numzgrid),stat=stat)
+  if (stat.ne.0) write(*,*)'ERROR: could not allocate volume'
+  allocate(areaeast(0:numxgrid-1,0:numygrid-1,numzgrid),stat=stat)
+  if (stat.ne.0) write(*,*)'ERROR: could not allocate areaeast'
+  allocate(areanorth(0:numxgrid-1,0:numygrid-1,numzgrid),stat=stat)
+  if (stat.ne.0) write(*,*)'ERROR: could not allocate areanorth'
   return
 
 
-999   write(*,*) ' #### FLEXPART MODEL ERROR! FILE "OUTGRID"    #### '
+999 write(*,*) ' #### FLEXPART MODEL ERROR! FILE "OUTGRID"    #### '
   write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
-  write(*,*) ' #### xxx/flexpart/options                    #### '
+  write(*,'(a)') path(1)(1:length(1))
   stop
 
+1000 write(*,*) ' #### FLEXPART MODEL ERROR! FILE "OUTGRID"    #### '
+  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
+  write(*,'(a)') path(2)(1:length(2))
+  stop
+
 end subroutine readoutgrid