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


Ignore:
Timestamp:
Sep 11, 2018, 6:06:31 PM (6 years ago)
Author:
pesei <petra seibert at univie ac at>
Branches:
univie
Children:
93786a1
Parents:
34f1452
Message:

modify most input read subroutines

changed some variable names (mostly for I-N reasons)
includes two names appearing also in timemanager, com_mod
corrected a few mistakes
simplified some parts of code
changed options/RELEASES which is in nml fmt correspondingly

File:
1 edited

Legend:

Unmodified
Added
Removed

  • src/readoutgrid.f90

    r8a65cb0 rd7935de  
    2929  !                                                                            *
    3030  !     4 June 1996                                                            *
    31   !     HSO, 1 July 2014
    32   !     Added optional namelist input
     31  !     HSO, 1 July 2014: Add optional namelist input                          *
     32  !     PS, 6/2015-9/2018: read regular input with free format                 *
     33  !       simplify code and rename some variables                              *
    3334  !                                                                            *
    3435  !*****************************************************************************
    … …  
    5152  implicit none
    5253
    53   integer :: i,j,stat
    54   real :: outhelp,xr,xr1,yr,yr1
     54  integer :: i,kz,istat
     55  real :: xr,xr1,yr,yr1
    5556  real,parameter :: eps=1.e-4
    5657
    5758  ! namelist variables
    5859  integer, parameter :: maxoutlev=500
    59   integer :: readerror
    60   real,allocatable, dimension (:) :: outheights
     60  integer :: ios
     61  real,allocatable, dimension (:) :: outheights,outaux
     62  logical :: lnml
    6163
    6264  ! declare namelist
    63   namelist /outgrid/ &
     65  namelist /nml_outgrid/ &
    6466    outlon0,outlat0, &
    6567    numxgrid,numygrid, &
    6668    dxout,dyout, &
    67     outheights
    68 
    69   ! allocate large array for reading input
    70   allocate(outheights(maxoutlev),stat=stat)
    71   if (stat.ne.0) write(*,*)'ERROR: could not allocate outheights'
    72 
    73   ! helps identifying failed namelist input
    74   dxout=-1.0
    75   outheights=-1.0
     69    outheight
     70
     71! allocate outheights for nml read with max dimension
     72  allocate(outheights(maxoutlev),outaux(maxoutlev),stat=istat)
     73  if (istat .ne. 0) write(*,*)'ERROR: could not allocate outheights'
    7674
    7775  ! Open the OUTGRID file and read output grid specifications
    7876  !**********************************************************
    7977
    80   open(unitoutgrid,file=path(1)(1:length(1))//'OUTGRID',status='old',form='formatted',err=999)
     78  outheight(:) = -999. ! initialise for later finding #valid levels
     79  open(unitoutgrid,file=trim(path(1))//'OUTGRID',status='old',&
     80    form='formatted',err=999)
    8181
    8282  ! try namelist input
    83   read(unitoutgrid,outgrid,iostat=readerror)
     83  read(unitoutgrid,nml_outgrid,iostat=ios)
    8484  close(unitoutgrid)
    8585
    86   if ((dxout.le.0).or.(readerror.ne.0)) then
    87 
    88     readerror=1
    89 
    90     open(unitoutgrid,file=path(1)(1:length(1))//'OUTGRID',status='old',err=999)
    91 
     86  if (ios .eq. 0) then ! namelist works
     87
     88    lnml = .true.
     89
     90  else ! read as regular text
     91
     92    lnml = .false.
     93
     94    open(unitoutgrid,file=trim(path(1))//'OUTGRID',status='old',err=999)
    9295    call skplin(5,unitoutgrid)
    9396
    94     ! 1. Read horizontal grid specifications
     97   ! Read horizontal grid specifications
    9598    !****************************************
    9699
    97100    call skplin(3,unitoutgrid)
    98     read(unitoutgrid,'(4x,f11.4)') outlon0
    99     call skplin(3,unitoutgrid)
    100     read(unitoutgrid,'(4x,f11.4)') outlat0
    101     call skplin(3,unitoutgrid)
    102     read(unitoutgrid,'(4x,i5)') numxgrid
    103     call skplin(3,unitoutgrid)
    104     read(unitoutgrid,'(4x,i5)') numygrid
    105     call skplin(3,unitoutgrid)
    106     read(unitoutgrid,'(4x,f12.5)') dxout
    107     call skplin(3,unitoutgrid)
    108     read(unitoutgrid,'(4x,f12.5)') dyout
    109 
    110   endif
     101    read(unitoutgrid,*) outlon0
     102    call skplin(3,unitoutgrid)
     103    read(unitoutgrid,*) outlat0
     104    call skplin(3,unitoutgrid)
     105    read(unitoutgrid,*) numxgrid
     106    call skplin(3,unitoutgrid)
     107    read(unitoutgrid,*) numygrid
     108    call skplin(3,unitoutgrid)
     109    read(unitoutgrid,*) dxout
     110    call skplin(3,unitoutgrid)
     111    read(unitoutgrid,*) dyout
     112
     113  endif ! read OUTGRID file
    111114
    112115  ! Check validity of output grid (shall be within model domain)
    … …  
    124127    write(*,*) ' #### GRID IS OUTSIDE MODEL DOMAIN. CHANGE    ####'
    125128    write(*,*) ' #### FILE OUTGRID IN DIRECTORY               ####'
    126     write(*,'(a)') path(1)(1:length(1))
     129    write(*,'(a)') trim(path(1))
    127130    stop
    128131  endif
    129132
    130   ! 2. Count Vertical levels of output grid
    131   !****************************************
    132 
    133   if (readerror.ne.0) then
    134     j=0
    135 100 j=j+1
    136     do i=1,3
    137       read(unitoutgrid,*,end=99)
    138     end do
    139     read(unitoutgrid,'(4x,f7.1)',end=99) outhelp
    140     if (outhelp.eq.0.) goto 99
    141     goto 100
    142 99  numzgrid=j-1
    143   else
    144     do i=1,maxoutlev
    145       if (outheights(i).lt.0) exit
    146     end do
    147     numzgrid=i-1
    148   end if
    149 
    150   allocate(outheight(numzgrid),stat=stat)
    151   if (stat.ne.0) write(*,*)'ERROR: could not allocate outheight'
    152   allocate(outheighthalf(numzgrid),stat=stat)
    153   if (stat.ne.0) write(*,*)'ERROR: could not allocate outheighthalf'
    154 
    155   ! 2. Vertical levels of output grid
    156   !**********************************
    157 
    158   if (readerror.ne.0) then
    159 
    160     rewind(unitoutgrid)
    161     call skplin(29,unitoutgrid)
    162 
    163     do j=1,numzgrid
    164       do i=1,3
    165         read(unitoutgrid,*)
    166       end do
    167       read(unitoutgrid,'(4x,f7.1)') outhelp
    168       outheight(j)=outhelp
    169       outheights(j)=outhelp
    170     end do
    171     close(unitoutgrid)
    172 
    173   else
    174 
    175     do j=1,numzgrid
    176       outheight(j)=outheights(j)
    177     end do
    178 
     133! Read (if .not. lmnl) and count vertical levels of output grid
     134!**************************************************************
     135
     136  do kz = 1,maxoutlev
     137    if (lnml) then ! we have read them already
     138      if (outheight(kz) .lt. 0.) exit ! 1st nondefined level
     139    else
     140      call skplin(3,unitoutgrid)
     141      read(unitoutgrid,*,end=10) outheight(kz)
     142    endif
     143  end do
     14410 continue 
     145
     146  numzgrid = kz - 1 ! number of outgrid levels
     147
     148! allocate the required length only, shuffle data 
     149  outaux = outheight ! shuffle
     150
     151  deallocate(outheights)
     152  allocate(outheight(numzgrid),outheighthalf(numzgrid),stat=istat)
     153  if (istat .ne. 0) then
     154    write(*,*) 'ERROR: could not allocate outheight and outheighthalf'
     155    stop 'readoutgrid error'
    179156  endif
    180157
    181   ! write outgrid file in namelist format to output directory if requested
    182   if (nmlout.and.lroot) then
    183     ! reallocate outheights with actually required dimension for namelist writing
    184     deallocate(outheights)
    185     allocate(outheights(numzgrid),stat=stat)
    186     if (stat.ne.0) write(*,*)'ERROR: could not allocate outheights'
    187 
    188     do j=1,numzgrid
    189       outheights(j)=outheight(j)
    190     end do
    191 
    192     open(unitoutgrid,file=path(2)(1:length(2))//'OUTGRID.namelist',err=1000)
    193     write(unitoutgrid,nml=outgrid)
     158  outheight=outaux(1:numzgrid) ! shuffle back
     159  deallocate (outaux)
     160
     161! write outgrid file in namelist format to output directory if requested
     162  if (nmlout) then
     163    open(unitoutgrid,file=trim(path(2))//'OUTGRID.namelist',err=1000)
     164    write(unitoutgrid,nml=nml_outgrid)
    194165    close(unitoutgrid)
    195166  endif
    … …  
    198169  !***************************************************************
    199170
    200   do j=2,numzgrid
    201     if (outheight(j).le.outheight(j-1)) then
    202     write(*,*) ' #### FLEXPART MODEL ERROR! YOUR SPECIFICATION#### '
    203     write(*,*) ' #### OF OUTPUT LEVELS IS CORRUPT AT LEVEL    #### '
    204     write(*,*) ' #### ',j,'                              #### '
    205     write(*,*) ' #### PLEASE MAKE CHANGES IN FILE OUTGRID.    #### '
    206     endif
     171  do kz=2,numzgrid
     172    if (outheight(kz) .le. outheight(kz-1)) goto 998
    207173  end do
    208174
    … …  
    210176  !*****************************************************************
    211177
    212   outheighthalf(1)=outheight(1)/2.
    213   do j=2,numzgrid
    214     outheighthalf(j)=(outheight(j-1)+outheight(j))/2.
     178  outheighthalf(1) = 0.5*outheight(1)
     179  do kz = 2,numzgrid
     180    outheighthalf(kz) = 0.5*(outheight(kz-1)+outheight(kz))
    215181  end do
    216182
    … …  
    218184  youtshift=ylat0-outlat0
    219185
    220   allocate(oroout(0:numxgrid-1,0:numygrid-1),stat=stat)
    221   if (stat.ne.0) write(*,*)'ERROR: could not allocate oroout'
    222   allocate(area(0:numxgrid-1,0:numygrid-1),stat=stat)
    223   if (stat.ne.0) write(*,*)'ERROR: could not allocate area'
    224   allocate(volume(0:numxgrid-1,0:numygrid-1,numzgrid),stat=stat)
    225   if (stat.ne.0) write(*,*)'ERROR: could not allocate volume'
    226   allocate(areaeast(0:numxgrid-1,0:numygrid-1,numzgrid),stat=stat)
    227   if (stat.ne.0) write(*,*)'ERROR: could not allocate areaeast'
    228   allocate(areanorth(0:numxgrid-1,0:numygrid-1,numzgrid),stat=stat)
    229   if (stat.ne.0) write(*,*)'ERROR: could not allocate areanorth'
     186  allocate(oroout(0:numxgrid-1,0:numygrid-1),stat=istat)
     187  if (istat .ne. 0) write(*,*)'ERROR: could not allocate oroout'
     188  allocate(area(0:numxgrid-1,0:numygrid-1),stat=istat)
     189  if (istat .ne. 0) write(*,*)'ERROR: could not allocate area'
     190  allocate(volume(0:numxgrid-1,0:numygrid-1,numzgrid),stat=istat)
     191  if (istat .ne. 0) write(*,*)'ERROR: could not allocate volume'
     192  allocate(areaeast(0:numxgrid-1,0:numygrid-1,numzgrid),stat=istat)
     193  if (istat .ne. 0) write(*,*)'ERROR: could not allocate areaeast'
     194  allocate(areanorth(0:numxgrid-1,0:numygrid-1,numzgrid),stat=istat)
     195  if (istat .ne. 0) write(*,*)'ERROR: could not allocate areanorth'
     196 
    230197  return
    231198
     199998 continue
     200  write(*,*) ' #### FLEXPART MODEL ERROR! YOUR SPECIFICATION#### '
     201  write(*,*) ' #### OF OUTPUT LEVELS NOT INCREASING AT LEVEL#### '
     202  write(*,*) ' #### ',kz,'                                  #### '
     203  write(*,*) ' #### PLEASE MAKE CHANGES IN FILE OUTGRID.    #### '
     204  STOP 'readoutgrid error'
    232205
    233206999 write(*,*) ' #### FLEXPART MODEL ERROR! FILE "OUTGRID"    #### '
    234207  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
    235   write(*,'(a)') path(1)(1:length(1))
     208  write(*,'(a)') trim(path(1))
    236209  stop
    237210
    2382111000 write(*,*) ' #### FLEXPART MODEL ERROR! FILE "OUTGRID"    #### '
    239212  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
    240   write(*,'(a)') path(2)(1:length(2))
     213  write(*,'(a)') trim(path(2))
    241214  stop
    242215
Note: See TracChangeset for help on using the changeset viewer.
hosted by ZAMG