EnKF/EnKF-MPI-TOPAZ/mod_analysisfields.F90

!KAL -- this module allows us to fine-tune the fields
!KAL -- we wish to include in the analysis. The new
!KAL -- layout of the EnKF makes it possible to specify fields
!KAL -- to analyze at run-time rather than at compile-time
!KAL -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
!KAL --
!KAL -- Module variables:
!KAL --    numfields   - total number of fields to process
!KAL --    fieldnames  - the names of the fields we wish to analyze
!KAL --    fieldlevel  - the levels of the associated fields
!KAL --    fieldtype   - in which file the field can be found: 
!ckb --                  1: ice, 2: ocean, 3: ice parameter, 
!ckb --                  4: ocean parameter
!KAL --
!KAL -- Ex: If we only want to assimilate temperatures in layer
!KAL --     one and two, numfields, fieldnames and fieldlevel 
!KAL --     would look like:
!KAL --
!KAL --     numfields=2                                 
!KAL --     fieldnames (1)='temp', fieldnames (2)='temp'
!KAL --     fieldlevel (1)=     1, fieldlevel (2)=2     (???)
!ckb --     fieldtype  (1)=     2, fieldtype  (2)=2
!KAL -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
!KAL -- The file "analysisfields.in" specifies the fields to 
!KAL -- inlude in the analysis. Format of one line is fieldname
!ckb -- first layer, last layer and fieldtype. For example:
!KAL --
!KAL -- fieldname       1 31  1
!KAL -- 12345678901234567890123456789012345678901234567890
!KAL --
!KAL -- Fortran format for one line is '(a14,3i3)'
!KAL --
!KAL -- Example: to specify that we want temperature and salinity 
!ckb --          in layers 1-31 (ocean variables, type 2) to be
!ckb --           updated, as well as ice concentration (layer 0,
!ckb --           type 1), and the atmosphere-ice-drag coefficient,
!ckb --           specify:
!ckb --
!ckb -- a_i_htc1        0  0  1
!ckb -- v_i_htc1        0  0  1
!ckb -- tempt_il3_htc5  0  0  1
!ckb -- ub              1 31  2
!ckb -- vb              1 31  2
!ckb -- cai             0  0  3
!KAL -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --

! [FM,CKB] Changed to allow for column "parameter"





module mod_analysisfields

character(len=*), parameter :: infile='analysisfields.in'
integer,save :: numfields
character(len=14), dimension(:), save, allocatable:: fieldnames
integer          , dimension(:), save, allocatable:: fieldlevel 
integer          , dimension(:), save, allocatable:: fieldtype

contains

   integer function get_nrfields()

#if defined (QMPI)
   use qmpi
#else
   use qmpi_fake
#endif

   implicit none
   integer :: ios,first,last,type
   logical :: ex
   character(len=14) :: char14

   inquire(exist=ex,file=infile)
   if (.not. ex) then
      if (master) print *,'Could not find '//infile
      call stop_mpi()
   end if

   open(10,status='old',form='formatted',file=infile)
   ios=0
   get_nrfields=0
   do while (ios==0)
      read(10,100,iostat=ios) char14,first,last,type
      if (ios==0) get_nrfields=get_nrfields+last-first+1
   end do
   close(10)
   100 format (a14,3i3)
   end function



   subroutine get_analysisfields()

#if defined (QMPI)
   use qmpi
#else
   use qmpi_fake
#endif

   implicit none
   integer :: first,last,type,k,nfld,ios
   logical :: ex
   character(len=14) :: char14

   numfields=get_nrfields()
   if (master) print *,'numfields is ',numfields
   if (numfields<=0 .or.numfields > 18000) then ! FM I Changed 600 to 18000
      if (master) print *,'(get_analysisfields) numfields is higher than max allowed setting or = 0'
      call stop_mpi()
   end if
   allocate(fieldnames(numfields))
   allocate(fieldlevel(numfields))
   allocate(fieldtype(numfields))


   inquire(exist=ex,file=infile)
   if (.not. ex) then
      if (master) print *,'Could not find '//infile
      call stop_mpi()
   end if

   open(10,status='old',form='formatted',file=infile)
   ios=0
   nfld=0
   do while (ios==0)
      read(10,100,iostat=ios) char14,first,last,type
      if (ios==0) then
         do k=first,last
            fieldnames (nfld+k-first+1)=char14
            fieldlevel (nfld+k-first+1)=k
            fieldtype  (nfld+k-first+1)=type
         end do
         nfld=nfld+last-first+1
      end if
   end do
   close(10)
   100 format (a14,3i3)

   if (nfld/=numfields) then
      if (master) print *,'An error occured when reading '//infile
      call stop_mpi()
   end if

   ! List fields used in analysis
   print *, "(mod_analysisfields) Fields used in analysis:"
   print *, "(mod_analysisfields) --- removed to reduce output ---"
   !do k=1,numfields
   !   if (master) print *,fieldnames(k),fieldlevel(k),fieldtype(k)
   !end do

   end subroutine
end module mod_analysisfields