  !> @brief Functions to parse command line arguments
  !>
  !> Module to read command line arguments to a program
  !> We assume they are of the form name=value. There must
  !> NOT be spaces around the '=' sign!
  !> Value can be extracted as a string, a long-integer
  !> or a double-precision real.
  !> Argument names are limited to 20 chars, and values
  !> to 30 chars as read.
  !>
  !> Note that the only functions you should call from outside are
  !> parse_args, get_arg and get_arg_value
  !> A complete example code is:
  ! @snippet command_line_snippet.f90 Cmd eg
  !> @include command_line_snippet.f90

  !> @author H Ratcliffe

MODULE command_line

  USE kinds

  IMPLICIT NONE

  !> Length of character value string
  INTEGER, PARAMETER :: len = 30

  !> Type containing a key-value pair
  ! Init. to default values
  TYPE cmd_arg
    CHARACTER(LEN=20) :: name = "NULL"
    CHARACTER(LEN=len) :: value = ""
  END TYPE

  ! Interface dispatches to the correct actual function
  ! based on arguments in call
  ! This means we can convert to a requested type
  ! We call like e.g.: get_arg(10, var) or
  ! get_arg("name", var) and we will try and
  ! interpret the value as whatever type var is
  ! We can supply optional arg to those
  ! to capture whether or not the name exists
  ! Note that if a name was present, but the value could
  ! not be parsed, exists will be TRUE
  ! but correct will be FALSE
  !> @brief Read arguments by name or number
  !>
  !> All of the members of this interface take either a 'name' or a 'num'
  !>  parameter, detailling which argument to look up and
  !> the type of the 'val' argument dictates the type
  !> to attempt parsing the argument as
  !
  !> @param name Name to look up (if used)
  !> @param num Arg. number to look up (if used)
  !> @param val Value to read into
  !> @param exists Whether the name was found
  !> @return True if the name is found and parsed, False otherwise
  INTERFACE get_arg
    MODULE PROCEDURE get_arg_num_int, get_arg_name_int
    MODULE PROCEDURE get_arg_num_dbl, get_arg_name_dbl
    MODULE PROCEDURE get_arg_num_str, get_arg_name_str
  END INTERFACE

  ! Module level variables - these are accessible
  ! to all functions in the module
  ! but we make them private to within the module only
  !> The argument list
  TYPE(cmd_arg), DIMENSION(:), ALLOCATABLE :: all_args
  !> The number of arguments
  INTEGER :: num_args = 0
  PRIVATE :: all_args, num_args

  CONTAINS

  !> @brief Parse out command line args
  ! This function can be called multiple times
  ! and will freshly parse ALL arguments each time
  ! We assume these are entered as 'name=value' and
  ! if there is no '=' sign, then value is set to a sentinel
  ! For most flexibility, we assume all values are reals, and
  ! thus we use 'HUGE(1.0)' as the no-value sentinel. This number
  ! is greater than any other real value
  SUBROUTINE parse_args()
    !Grab all command line args and store

    ! Strictly we can't be sure 50 chars is enough
    ! but for command-line args it's enough if we're sensible
    ! We wont overflow, but our strings may get truncated
    CHARACTER(LEN=51) :: arg
    INTEGER :: i, indx

    num_args = COMMAND_ARGUMENT_COUNT()
    IF(num_args > 0) THEN

      ! If this is not the first call, all_args may be already allocated
      ! Deallocate if needed, and allocate to correct size
      IF(ALLOCATED(all_args)) DEALLOCATE(all_args)
      ALLOCATE(all_args(num_args))

      ! Loop over all arguments
      DO i = 1, num_args
        CALL GET_COMMAND_ARGUMENT(i, arg)
        ! Location of the '=' sign
        ! If not found, return value is 0
        indx = INDEX(arg, '=')
        IF(indx > 1) THEN
          ! All characters up to '=', not including it
          ! but with any leading spaces removed
          all_args(i)%name = ADJUSTL(arg(1:indx-1))
          ! All characters after '='
          all_args(i)%value= ADJUSTL(arg(indx+1:))
       ELSE
          all_args(i)%name = TRIM(ADJUSTL(arg))
          ! Value already has a default value, so leave it alone
        END IF
      END DO
    ENDIF

  END SUBROUTINE parse_args

!------------------------------------------------------------------
  ! The next functions let you access the parsed args
  ! You need to have called parse_args first or
  ! there wont be any args!
  ! You don't need to call these explicitly, you should
  ! just use get_arg

  !> @brief Read by number for double precision values
  !> @param num Argument number to read
  !> @param val Value to read into
  !> @param exists Whether the name was found
  !> @return True if the name is found and parsed, False otherwise
  FUNCTION get_arg_num_dbl(num, val, exists)

    LOGICAL :: get_arg_num_dbl
    INTEGER, INTENT(IN) :: num
    REAL(KIND=REAL64), INTENT(OUT) :: val
    LOGICAL, INTENT(OUT), OPTIONAL :: exists
    LOGICAL :: found
    INTEGER :: ierr

    found = .FALSE.
    ! Check requested number is in range
    IF(num <= num_args .AND. num > 0) THEN
      ! READ it from string into value
      ! We don't need to specify the format in general
      READ(all_args(num)%value, *, IOSTAT=ierr) val
      found = .TRUE.
    END IF

    IF(PRESENT(exists)) THEN
      exists = found
    END IF

    ! Return value is whether value is found and correctly parsed
    get_arg_num_dbl = (found .AND. (ierr == 0))

  END FUNCTION get_arg_num_dbl

  !> @brief Read by name for double precision values
  !> @param name Argument name to look up
  !> @param val Value to read into
  !> @param exists Whether the name was found
  !> @return True if the name is found and parsed, False otherwise
  FUNCTION get_arg_name_dbl(name, val, exists)

    LOGICAL :: get_arg_name_dbl
    CHARACTER(LEN=*), INTENT(IN) :: name
    REAL(KIND=REAL64), INTENT(OUT) :: val
    INTEGER :: i
    LOGICAL, INTENT(OUT), OPTIONAL :: exists
    LOGICAL :: found
    INTEGER :: ierr

    found = .FALSE.
    ! Our cmd_arg type is already initialised to the sentinel
    DO i = 1, num_args
      IF(all_args(i)%name == TRIM(ADJUSTL(name))) THEN
        found = .TRUE.
        READ(all_args(i)%value, *, IOSTAT=ierr) val
        EXIT
      END IF
    END DO

    IF(PRESENT(exists)) THEN
      exists = found
    END IF

    ! Return value is whether value is found and correctly parsed
    get_arg_name_dbl = (found .AND. (ierr == 0))

  END FUNCTION get_arg_name_dbl

  !> @brief Read by number for long integer values
  !> @param num Argument number to read
  !> @param val Value to read into
  !> @param exists Whether the name was found
  !> @return True if the name is found and parsed, False otherwise
  FUNCTION get_arg_num_int(num, val, exists)

    LOGICAL :: get_arg_num_int
    INTEGER, INTENT(IN) :: num
    INTEGER(KIND=INT32), INTENT(OUT) :: val
    LOGICAL, INTENT(OUT), OPTIONAL :: exists
    LOGICAL :: found
    INTEGER :: ierr

    found = .FALSE.
    ! Check requested number is in range
    IF(num <= num_args .AND. num > 0) THEN
      ! READ it from string into value
      ! We don't need to specify the format in general
      READ(all_args(num)%value, *, IOSTAT=ierr) val
      found = .TRUE.
    END IF

    IF(PRESENT(exists)) THEN
      exists = found
    END IF

    ! Return value is whether value is found and correctly parsed
    get_arg_num_int = (found .AND. (ierr == 0))

  END FUNCTION get_arg_num_int

  !> @brief Read by name for long integer values
  !> @param name Argument name to look up
  !> @param val Value to read into
  !> @param exists Whether the name was found
  !> @return True if the name is found and parsed, False otherwise
  FUNCTION get_arg_name_int(name, val, exists)

    LOGICAL :: get_arg_name_int
    CHARACTER(LEN=*), INTENT(IN) :: name
    INTEGER(KIND=INT32), INTENT(OUT) :: val
    INTEGER :: i
    LOGICAL, INTENT(OUT), OPTIONAL :: exists
    LOGICAL :: found
    INTEGER :: ierr

    found = .FALSE.
    ! Our cmd_arg type is already initialised to the sentinel
    DO i = 1, num_args
      IF(all_args(i)%name == TRIM(ADJUSTL(name))) THEN
        found = .TRUE.
        READ(all_args(i)%value, *, IOSTAT=ierr) val
        EXIT
      END IF
    END DO

    IF(PRESENT(exists)) THEN
      exists = found
    END IF

    ! Return value is whether value is found and correctly parsed
    get_arg_name_int = (found .AND. (ierr == 0))

  END FUNCTION get_arg_name_int

  !> @brief Read by number for string/character values
  !> @param num Argument number to read
  !> @param val Value to read into
  !> @param exists Whether the name was found - this is already contained in
  !> the return value, but is given for consistency with the other members
  !> @return True if the name is found and parsed, False otherwise
  FUNCTION get_arg_num_str(num, val, exists)

    LOGICAL :: get_arg_num_str
    INTEGER, INTENT(IN) :: num
    CHARACTER(LEN=*), INTENT(OUT) :: val
    LOGICAL, INTENT(OUT), OPTIONAL :: exists
    LOGICAL :: found

    found = .FALSE.
    ! Check requested number is in range
    IF(num <= num_args .AND. num > 0) THEN
      ! READ it from string into value
      ! We don't need to specify the format in general
      val = all_args(num)%value
      found = .TRUE.
    END IF

    IF(PRESENT(exists)) THEN
      exists = found
    END IF

    ! Return value is whether value is found and correctly parsed
    get_arg_num_str = found

  END FUNCTION get_arg_num_str

  !> @brief Read by name for string values
  !> @param name Argument name to look up
  !> @param val Value to read into
  !> @param exists Whether the name was found - this is already contained in
  !> the return value, but is given for consistency with the other members
  !> @return True if the name is found and parsed, False otherwise
  FUNCTION get_arg_name_str(name, val, exists)

    LOGICAL :: get_arg_name_str
    CHARACTER(LEN=*), INTENT(IN) :: name
    CHARACTER(LEN=*), INTENT(OUT) :: val
    INTEGER :: i
    LOGICAL, INTENT(OUT), OPTIONAL :: exists
    LOGICAL :: found

    found = .FALSE.
    ! Our cmd_arg type is already initialised to the sentinel
    DO i = 1, num_args
      IF(all_args(i)%name == TRIM(ADJUSTL(name))) THEN
        found = .TRUE.
        val = all_args(i)%value
        EXIT
      END IF
    END DO

    IF(PRESENT(exists)) THEN
      exists = found
    END IF

    ! Return value is whether value is found and correctly parsed
    get_arg_name_str = found

  END FUNCTION get_arg_name_str

!--------------------------------------------------------------------

  ! These lets you get just the string value from an argument by name
  ! you still need to have called parse_args first!
  !> @brief Lookup an argument by name and return the value as a string
  !> @param name Argument name to look up
  !> @param exists Whether the name was found
  !> @return The string value associated with the given name
  FUNCTION get_arg_value(name, exists)

    CHARACTER(LEN=len) :: get_arg_value
    CHARACTER(LEN=*), INTENT(IN) :: name
    LOGICAL, INTENT(OUT), OPTIONAL :: exists
    TYPE(cmd_arg) :: tmp
    INTEGER :: i
    LOGICAL :: found

    ! Initialise to the default value
    ! Use a temporary to get the default values
    ! This makes sure we match the expected sentinel
    get_arg_value = tmp%value
    found = .FALSE.

    DO i = 1, num_args
      IF(all_args(i)%name .EQ. TRIM(ADJUSTL(name))) THEN
        get_arg_value = all_args(i)%value
        found = .TRUE.
        EXIT
      END IF
    END DO

    IF(PRESENT(exists)) THEN
      exists = found
    END IF

  END FUNCTION get_arg_value

END MODULE command_line


text:
! Copyright (c) 1999, NVIDIA CORPORATION.  All rights reserved.
!
! Licensed under the Apache License, Version 2.0 (the "License");
! you may not use this file except in compliance with the License.
! You may obtain a copy of the License at
!
!     http://www.apache.org/licenses/LICENSE-2.0
!
! Unless required by applicable law or agreed to in writing, software
! distributed under the License is distributed on an "AS IS" BASIS,
! WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
! See the License for the specific language governing permissions and
! limitations under the License.
!
! Test noncontiguous pointers.
! Bug when an outer-scope pointer is referenced in an internal procedure
!
integer, parameter:: N=10
integer, target :: a(10,N)
integer, pointer :: p(:)
integer :: expect(N) = (/1,2,3,4,5,6,7,8,9,10/)
integer :: result(N)

p=>a(2,:)
call foo

!print *,  a(2,:)
result = a(2,1:N)
call check(result, expect, N)

contains
    subroutine foo
    do i = 1, N
        p(i) = i   !p's SDSCNS1 flag must be set in the backend
    enddo
    endsubroutine
end

max_stars_repo_name: kammerdienerb/flang

max_issues_repo_name: kammerdienerb/flang

max_forks_repo_name: kammerdienerb/flang


text:
      double precision function xstar(eps1,eps2, vec1,vec2, ndeg, elpty,
     &       ilinit)
      implicit double precision (a-h, o-z)
c     calculating nstar=deg*cos(eps r1)*cos(eps rN)
c     written by alexei ankudinov 08.13.96
c     calculate the plane wave approximation for central atom
c     vec1 - direction to the first atom in path
c     vec2 - direction to the last atom in path
c     ndeg - may be not equal to 'deg' in diff. paths
c             subroutines

      include '../HEADERS/const.h'
      include '../HEADERS/dim.h'
c      include 'pdata.h'
      integer ilinit

      dimension eps1(3), eps2(3), vec1(3), vec2(3)

      lfin = ilinit
      x = xxcos(vec1, vec2)
      iav = 1
      y = xxcos(eps1,vec1)
      z = xxcos(eps1,vec2)
      xtemp = ystar(lfin, x, y, z, iav)
      if (elpty .ne. 0.0) then
        y = xxcos(eps2,vec1)
        z = xxcos(eps2,vec2)
        xtemp = xtemp + elpty**2 * ystar(lfin, x, y, z, iav)
      endif
      xstar = ndeg * xtemp /(1+elpty**2)

      return
      end

      double precision function xxcos (veca, vecb)
      implicit double precision (a-h, o-z)
      dimension veca(3), vecb(3)

      x1 = 0
      do 23 j = 1,3
         x1 = x1 + veca(j) * vecb(j)
   23 continue
      xnorma = 0
      xnormb = 0
      do 24 j = 1,3
         xnorma = xnorma + veca(j)**2
         xnormb = xnormb + vecb(j)**2
   24 continue
      xxcos = x1/sqrt(xnorma*xnormb)
      return
      end

      double precision function ystar (lfin, x , y, z, iav)
      implicit double precision (a-h, o-z)
c
      dimension pln (0:4,4)
      data pln /0.0  , 1.0, 0.0 , 0.0, 0.0,
     2         -0.5  , 0.0, 1.5 , 0.0, 0.0,
     3          0.0  ,-1.5, 0.0 , 2.5, 0.0,
     4          0.375, 0.0,-3.75, 0.0, 4.375/

      pln0  = pln(0,lfin)
      do 40 i = 1, lfin
         pln0  = pln0  + pln(i, lfin) * x**i
  40  continue
      if (iav.eq.0) then
         ystar = pln0/(2*lfin+1)
      else
         pln1  = pln(1,lfin)
         do 50 i = 2, lfin
            pln1  = pln1  + pln(i, lfin)*i*x**(i-1)
  50     continue
         pln2  = 2* pln(2,lfin)
         do 60 i = 3, lfin
            pln2  = pln2  + pln(i, lfin)*i*(i-1)*x**(i-2)
  60     continue
         ytemp = - lfin*pln0 + pln1*(x+y*z) - pln2*(y**2+z**2-2*x*y*z)
         ystar = ytemp * 3/lfin/(4*lfin**2-1)
      endif
      return
      end

ext: f


MODULE DFCOMM

  use types
  implicit none
!      This module allocated arrays for the opacity pre-tabulation calcualtions!
!
!.... maxd   = maximum number of depths in the atmosphere
!.... maxmeq = maximum number of molecular equations
!.... maxloc = maximum number of molecular components
!.... maxmol = maximum number of molecules
!.... maxmu  = maximum number of angles for the radiation field
!.... maxnu  = the maximum number of frequencies

! --------------------------------------------

  integer, parameter :: maxd = 360
  integer, parameter :: lenbuff = 3510000
  integer, parameter :: maxprof = 10000
!
!  for ODF's T and P gird

  real(kind=dp) :: profile(10000), b10000(10000)
  real(kind=dp) :: extab(10001), e1tab(2000)
  real(kind=dp) :: h0tab(2001), h1tab(2001), h2tab(2001), atab(281,1001)

  real(kind=dp) :: xnfpel(594),dopple(594),xnfdop(594) !6*99=594 -> 6 levels and 99 atoms
  real(kind=dp) :: buffer(lenbuff), continuum(lenbuff)

 ! for  Kuruzc  bin-sizes
  real(kind=dp) :: biga(95),bigb(125),bigc(109)
  real(kind=dp) :: lita(95),litb(190),litc(190),litd(190),lite(190)
  real(kind=dp) :: litf(190),litg(168)

  real(kind=dp) :: wlbeg,wlend, wlbeg0, resolu,ratio,ratiolg,wbegin, wend
  integer :: ixwlbeg, length, mlines
  integer :: nsizebig, nsizelit
  integer :: nufreq, nufreq3

! for dynamical binning
  integer :: nsubbin

!--- all allocatable arrays

  integer,  allocatable :: nbeg(:), nend(:)
  integer,  allocatable ::  nbigbeg(:), nbigend(:), nlitbeg(:), nlitend(:)
  integer,  allocatable :: nn(:), nsteps(:,:)
  integer(kind=2), allocatable :: insteps(:), iodfstep(:,:,:,:,:)
  integer(kind=2), allocatable :: iodfsendl(:,:,:,:), iodfrecvl(:,:,:,:)
  integer(kind=2), allocatable :: iodfsendb(:,:,:,:), iodfrecvb(:,:,:,:)

  integer, allocatable :: nvt(:)

  real(kind=dp), allocatable :: wavebig(:), wavelit(:)
  real(kind=dp), allocatable ::  emerge(:), txnxn(:), bstim(:), xnfh2(:)

  real(kind=dp), allocatable :: dopmax(:,:,:)

  real(kind=dp), allocatable :: continall(:,:), contscat(:,:)
  real(kind=dp), allocatable :: dopple4(: ,: ,:), xnfpel4(:,:,:)
  real(kind=dp), allocatable  :: vtprof(:,:),vtcenter(:)

! allocatable frequency grid
  real(kind=dp), allocatable  ::  inifreset(:)

! for dynamical binning
  real(kind=dp), allocatable :: binbound(:,:),sbwith(:,:)
!-----------------------------------------------------------------------------------
!clsa
!related to filters
      integer :: nnsize, Nlambda_phi, numssbin, nntot
      integer :: iphi, ic
      integer, allocatable :: nnsteps(:,:),indexx(:)
      real (kind=dp),allocatable :: wavbigtilda(:),wavlittilda(:)
      real (kind=dp), allocatable :: dlambda(:,:), dlamtilda(:,:), phi_hr(:,:)
      real (kind=dp), allocatable :: xlambda_phi(:),filter(:)
      real (kind=dp) :: phi_lp,phi_avg,dlam,lambda0,lambda,lambdatilda,lambdatildaus,totdlam
      real (kind=dp) :: xlambda_min, xlambda_max
!-------------------------------------------------------------------------------------------
      real (kind=dp), allocatable :: sb_weight_send(:,:), sb_weight_recv(:,:)
!
!-------------------------------------------------------------------------------------------
!-----------------------------------------------------------------------------------

CONTAINS

SUBROUTINE set_arrays(numsbin, numssbin,ntemp)

  implicit none
  include 'common.rhoxbl'
  include 'common.turbpr'
  include 'common.ifopbl'

  integer, intent(in) ::ntemp, numsbin
  integer, intent(in) :: numssbin
  integer :: ierr(10), i

 do i = 1, 10
  ierr(i)=-1
 end do

 if (numsbin .eq. 0) then
  nsizebig = 328
  nsizelit = 1212
 else
  nsizebig = numsbin-1
  nsizelit = 1
 endif

 allocate(iodfsendb(numssbin,nsizebig,nrhox,ntemp), stat = ierr(9))
 allocate(iodfrecvb(numssbin,nsizebig,nrhox,ntemp), stat = ierr(10))

 allocate(insteps(numssbin), stat=ierr(3))
 allocate(iodfstep(numssbin,nsizebig+nsizelit,nrhox,ntemp,numvt), stat = ierr(4) )

 allocate(nbeg(nsizebig+nsizelit), stat = ierr(5))
 allocate(nend(nsizebig+nsizelit), stat = ierr(6))
 allocate(nn(nsizebig+nsizelit), stat = ierr(7))
 allocate(nsteps(numssbin+2,nsizebig+nsizelit), stat= ierr(8))

 allocate(nbigbeg(nsizebig), stat = ierr(9))
 allocate(nbigend(nsizebig), stat = ierr(10))

 if (maxval(ierr) .gt. 0) then
   print*, 'first ten arrays could not be allocated'
   stop
 end if

!clsa
 nnsize=nsizebig+nsizelit
 if (ifilter .eq. 1) then
   allocate (indexx(lenbuff), stat=ierr(1))
   allocate (dlambda(nnsize,lenbuff), stat=ierr(2))
   allocate (dlamtilda(nnsize,lenbuff), stat= ierr(3))
   allocate (phi_hr(nnsize,lenbuff), stat=ierr(4))
   allocate (wavbigtilda(nsizebig+1), stat=ierr(5))
   allocate (wavlittilda(nsizelit+1), stat=ierr(6))
   allocate (nnsteps(numssbin+1,nnsize), stat=ierr(7))

   allocate(sb_weight_send(numssbin,nsizebig+nsizelit), stat = ierr(10))
   allocate(sb_weight_recv(numssbin,nsizebig+nsizelit), stat = ierr(2))

  if (maxval(ierr) .gt. 0) then
   print*, 'some of the arrays for filters could not be allocated'
   stop
 end if

 end if

!clsa

 allocate(nlitbeg(nsizelit), stat = ierr(1))
 allocate(nlitend(nsizelit), stat = ierr(2))

 allocate(wavebig(nsizebig+1), stat = ierr(3))
 allocate(wavelit(nsizelit+1), stat = ierr(4))

 allocate(binbound(nsizebig+nsizelit, numssbin),stat = ierr(5))
 allocate(sbwith(numssbin, nsizebig+nsizelit), stat = ierr(6))

 allocate(emerge(nrhox), stat=ierr(7))
 allocate(txnxn(nrhox), stat=ierr(8))
 allocate(bstim(nrhox), stat=ierr(9))

 allocate(xnfh2(nrhox), stat=ierr(10))

 if (maxval(ierr) .gt. 0) then
   print*, 'second ten arrays could not be allocated'
   stop
 end if

! allocate(tsave(nrhox))
 allocate(dopmax(nufreq3, 6, 99), stat=ierr(1))
 allocate(continall(nufreq, nrhox), stat=ierr(2))
 allocate(contscat(nufreq, nrhox), stat=ierr(3))
 allocate(dopple4(6,99, nrhox), stat=ierr(4))
 allocate(xnfpel4(6,99, nrhox), stat=ierr(5))

 allocate(vtprof(120,numvt), stat=ierr(6))
 allocate(vtcenter(numvt), stat=ierr(7))
 allocate(nvt(numvt),stat=ierr(8))

 allocate(iodfsendl(numssbin,nsizelit,nrhox,ntemp), stat = ierr(9))
 allocate(iodfrecvl(numssbin,nsizelit,nrhox,ntemp), stat = ierr(10))

 if (maxval(ierr) .gt. 0) then
   print*, 'second ten arrays could not be allocated'
   stop
 end if

 nbeg=0
 nend=0
 nn =0
 nsteps=0
 nbigbeg=0
 nbigend =0
 nlitbeg= 0
 nlitend= 0

 wavebig=0.0d0
 wavelit=0.0d0

 emerge=0.0d0
 txnxn = 0.0d0
 bstim = 0.0d0
 xnfh2 = 0.0d0
 dopmax = 0.0d0
 continall = 0.0d0
 contscat = 0.0d0
 dopple4 = 0.0d0
 xnfpel4 = 0.0d0

 iodfstep = 0
 iodfsendl = 0
 iodfrecvl = 0
 iodfsendb = 0
 iodfrecvb = 0

!lsa
 if (ifilter .eq. 1 ) then
  sb_weight_send=0.
  sb_weight_recv=0.
 end if

END SUBROUTINE set_arrays

SUBROUTINE close_arrays

  implicit none
  include 'common.ifopbl'

 deallocate(nbeg)
 deallocate(nend)
 deallocate(nn)
 deallocate(nsteps)

 deallocate(nbigbeg)
 deallocate(nbigend)
 deallocate(nlitbeg)
 deallocate(nlitend)

 deallocate(insteps)
 deallocate(iodfstep)
 deallocate(iodfsendb)
 deallocate(iodfrecvb)

 deallocate(iodfsendl)
 deallocate(iodfrecvl)

 deallocate(wavebig)
 deallocate(wavelit)

 deallocate(binbound)

 deallocate(emerge)
 deallocate(txnxn)
 deallocate(bstim)
 deallocate(xnfh2)
! deallocate(tsave)
 deallocate(dopmax)
 deallocate(continall)
 deallocate(contscat)
 deallocate(dopple4)
 deallocate(xnfpel4)

 deallocate(vtprof)
 deallocate(vtcenter)
 deallocate(nvt)

!lsa
 if (ifilter .eq. 1) then
   deallocate (indexx)
   deallocate (dlambda)
   deallocate (dlamtilda)
   deallocate (phi_hr)
   deallocate (wavbigtilda)
   deallocate (wavlittilda)
   deallocate (nnsteps)

   deallocate (sb_weight_send)
   deallocate (sb_weight_recv)
!lsa
 endif

END SUBROUTINE close_arrays

subroutine alloc_wave(numsbin, numssbin)

  implicit none
  integer, intent(in) :: numsbin, numssbin

  include 'common.rhoxbl'
  include 'common.turbpr'

 if (numsbin .eq. 0) then
  nsizebig = 328
  nsizelit = 1212

 else
  nsizebig = numsbin-1
  nsizelit = 1
 endif

 allocate(wavebig(nsizebig+1))
 allocate(wavelit(nsizelit+1))

end subroutine alloc_wave

subroutine close_wave

  implicit none
  deallocate(wavebig)
  deallocate(wavelit)

end subroutine close_wave

SUBROUTINE def_binsize

   implicit none
   include 'common.ifopbl'

   integer i

   if (ifkbin) then

     biga=(/8.97666,     9.2,         9.5,         9.71730,     9.81590, &
     &    10.10901,    10.3,        10.46451,    10.65,       10.88545, &
     &    11.2,        11.6,        11.95562,    12.3,        12.66565, &
     &    12.74676,    12.93504,    13.16033,    13.31723,    13.46700, &
     &    13.86253,    14.10462,    14.54782,    14.9,        15.3, &
     &    15.74874,    16.00346,    16.4,        16.8,        17.25088, &
     &    17.43258,    17.93660,    18.10114,    18.24438,    18.97796, &
     &    19.22422,    19.55098,    20.11987,    20.23592,    20.8, &
     &    21.3,        21.94783,    22.01983,    22.58070,    22.78377, &
     &    23.2,        23.65221,    24.3,        25.00744,    25.46614, &
     &    25.89239,    26.14316,    26.6,        27.1,        27.6, &
     &    28.09263,    28.65871,    29.3,        29.95849,    30.26260, &
     &    31.3,        32.3,        33.3,        34.3,        35.27931, &
     &    36.14100,    37.,         38.01600,    38.96642,    40., &
     &    41.,         41.87992,    42.5,        43.5,        44.73260, &
     &    45.5,        46.5,        47.5,        48.50178,    49.5, &
     &    50.42590,    50.87630,    52.,         54.,         56., &
     &    57.40614,    59.5,        61.5,        63.5,        65.53786, &
     &    67.,         69.,         71.,         72.24011,    74./)

    bigb=(/ 76.,         78.,         80.,         82.,         84., &
     &    86.,         88.,         90.,         91.17535,    94., &
     &    98.,         102.,       106.,        110.03056,   113., &
     &   116.,         120.,       123.92928,   128.,        132., &
     &   136.,         140.,       144.43391,   148.,        151.43485, &
     &   157.,         162.15072,  167.40282,   171.,        175., &
     &   180.,         184.,       188.,        193.,        197.46990, &
     &   202.,         207.13210,  210.,        215.,        220., &
     &   225.,         230.,       235.,        240.,        245., &
     &   251.12621,    255.,       260.,        265.,        270., &
     &   275.,         280.,       285.,        290.,        300., &
     &   310.,         320.,       330.,        340.,        350., &
     &   360.,         364.70183,  370.,        380.,        390., &
     &   400.,   410., 420., 430., 450., 460.,  470., 480.,  490., 500., &
     &   510.,   520., 530., 540., 550., 560.,  570., 580.,  590., 600., &
     &   610.,   620., 630., 640., 650., 660.,  670., 680.,  690., 700., &
     &   710.,   720., 730., 740., 750., 760.,  770., 780.,  790., 800., &
     &   810.,820.58271,830.,840., 850., 860.,  870., 880.,  890., 900., &
     &   910.,   920., 930., 940., 950., 960.,  970., 980.,  990.,1000./)
     bigc=(/1025.,1050.,1075.,1100.,1125.,1150.,1175.,1200.,1225.,1250., &
     & 1275.,1300.,1325.,1350.,1375.,1400.,1425.,1458.81670,1475.,1500., &
     & 1525.,1550.,1575.,1600.,1640.,1680.,1720.,1760.,1800.,1840., &
     & 1880.,1920.,1960.,2000.,2050.,2100.,2150.,2200.,2250.,2279.40330, &
     & 2300.,2350.,2400.,2450.,2500.,2550.,2600.,2650.,2700.,2750., &
     & 2800.,2850.,2900.,2950.,3000.,3050.,3100.,3150.,3200.,3282.34320, &
     & 3400.,3500.,3600.,3700.,3800.,3900.,4000.,4100.,4200.,4300., &
     & 4400.,4500.,4600.,4700.,4800.,4900.,5000.,5100.,5200.,5300., &
     & 5400.,5500.,5600.,5700.,5800.,5900.,6000.,6100.,6200.,6300., &
     & 6400.,6600.,6800.,7000.,7200.,7400.,7600.,7800.,8000.,8200., &
     & 8400.,8600.,8800.,9000.,9200.,9400.,9600.,9800.,10000./)

     lita=(/8.97666,     9.2,         9.5,         9.71730,     9.81590, &
     &    10.10901,    10.3,        10.46451,    10.65000,    10.88545, &
     &    11.2,        11.6,        11.95562,    12.3,        12.66565, &
     &    12.74676,    12.93504,    13.16033,    13.31723,    13.46700, &
     &    13.86253,    14.10462,    14.54782,    14.9,        15.3, &
     &    15.74874,    16.00346,    16.4,        16.8,        17.25088, &
     &    17.43258,    17.93660,    18.10114,    18.24438,    18.97796, &
     &    19.22422,    19.55098,    20.11987,    20.23592,    20.8, &
     &    21.3,        21.94783,    22.01983,    22.58070,    22.78377, &
     &    23.2,        23.65221,    24.3,        25.00744,    25.46614, &
     &    25.89239,    26.14316,    26.6,        27.1,        27.6, &
     &    28.09263,    28.65871,    29.3,        29.95849,    30.26260, &
     &    31.3,        32.3,        33.3,        34.3,        35.27931, &
     &    36.14100,    37.,         38.01600,    38.96642,    40., &
     &    41.,         41.87992,    42.5,        43.5,        44.73260, &
     &    45.5,        46.5,        47.5,        48.50178,    49.5, &
     &    50.42590,    50.87630,    51.5,        52.5,        53.5, &
     &    54.5,        55.5,        56.5,        57.40614,    58.5, &
     &    59.5,        60.5,        61.5,        62.5,        63.5/)
    litb=(/ 64.5, 65.53786, 66., 67., 68., 69., 70., 71., 72.24011,  73., &
     &  74., 75., 76., 77., 78., 79., 80., 81., 82., 83., &
     &  84., 85., 86., 87., 88., 89., 90., 91.17535, 92., 93., &
     &  94., 95., 96., 97., 98., 99., 100., 101., 102., 103., &
     & 104., 105., 106., 107., 108., 109., 110.03056, 111., 112., 113., &
     & 114., 115., 116., 117., 118., 119., 120., 121., 122., 123., &
     & 123.92928, 125., 126., 127., 128., 129., 130., 131., 132., 133., &
     & 134., 135., 136., 137., 138., 139., 140., 141., 142., 143., &
     & 144.,144.43391,145.,146.,147.,148.,149.,150.,151.,151.43485, &
     & 152., 153., 154., 155., 156., 157., 158., 159., 160., 161., &
     & 162.15072,163.,164.,165.,166.,167.,167.40282,168.,169.,170., &
     & 171., 172., 173., 174., 175., 176., 177., 178., 179., 180., &
     & 181., 182., 183., 184., 185., 186., 187., 188., 189., 190., &
     & 191.,192.,193.,194.,195.,196.,197.46990,197.81149,199.,200., &
     & 201.,202.,203.,204.,205.,206.,207.13210,207.61400,208.,209., &
     & 210., 211., 212., 213., 214., 215., 216., 217., 218., 219., &
     & 220., 221., 222., 223., 224., 225., 226., 227., 228., 229., &
     & 230., 231., 232., 233., 234., 235., 236., 237., 238., 239., &
     & 240., 241., 242., 243., 244., 245., 246., 247., 248., 249./)
    litc=(/250.,251.12621,251.51005,252.,253.,254.,255.,256.,257.,258., &
     & 259., 260., 261., 262., 263., 264., 265., 266., 267., 268., &
     & 269., 270., 271., 272., 273., 274., 275., 276., 277., 278., &
     & 279., 280., 281., 282., 283., 284., 285., 286., 287., 288., &
     & 289., 290., 292., 294., 296., 298., 300., 302., 304., 306., &
     & 308., 310., 312., 314., 316., 318., 320., 322., 324., 326., &
     & 328., 330., 332., 334., 336., 338., 340., 342., 344., 346., &
     & 348., 350., 352., 354., 356., 358., 360., 362., 364., 364.70183, &
     & 366., 368., 370., 372., 374., 376., 378., 380., 382., 384., &
     & 386., 388., 390., 392., 394., 396., 398., 400., 402., 404., &
     & 406., 408., 410., 412., 414., 416., 418., 420., 422., 424., &
     & 426., 428., 430., 432., 434., 436., 438., 440., 442., 444., &
     & 446., 448., 450., 452., 454., 456., 458., 460., 462., 464., &
     & 466., 468., 470., 472., 474., 476., 478., 480., 482., 484., &
     & 486., 488., 490., 492., 494., 496., 498., 500., 502., 504., &
     & 506., 508., 510., 512., 514., 516., 518., 520., 522., 524., &
     & 526., 528., 530., 532., 534., 536., 538., 540., 542., 544., &
     & 546., 548., 550., 552., 554., 556., 558., 560., 562., 564., &
     & 566., 568., 570., 572., 574., 576., 578., 580., 582., 584./)
    litd=(/586., 588., 590., 592., 594., 596., 598., 600., 602., 604., &
     & 606., 608., 610., 612., 614., 616., 618., 620., 622., 624., &
     & 626., 628., 630., 632., 634., 636., 638., 640., 642., 644., &
     & 646., 648., 650., 652., 654., 656., 658., 660., 662., 664., &
     & 666., 668., 670., 672., 674., 676., 678., 680., 682., 684., &
     & 686., 688., 690., 692., 694., 696., 698., 700., 702., 704., &
     & 706., 708., 710., 712., 714., 716., 718., 720., 722., 724., &
     & 726., 728., 730., 732., 734., 736., 738., 740., 742., 744., &
     & 746., 748., 750., 752., 754., 756., 758., 760., 762., 764., &
     & 766., 768., 770., 772., 774., 776., 778., 780., 782., 784., &
     & 786., 788., 790., 792., 794., 796., 798., 800., 802., 804., &
     & 806., 808., 810., 812., 814., 816., 818., 820.58271, 822., 824., &
     & 826., 828., 830., 832., 834., 836., 838., 840., 842., 844., &
     & 846., 848., 850., 852., 854., 856., 858., 860., 862., 864., &
     & 866., 868., 870., 872., 874., 876., 878., 880., 882., 884., &
     & 886., 888., 890., 892., 894., 896., 898., 900., 902., 904., &
     & 906., 908., 910., 912., 914., 916., 918., 920., 922., 924., &
     & 926., 928., 930., 932., 934., 936., 938., 940., 942., 944., &
     & 946., 948., 950., 952., 954., 956., 958., 960., 962., 964./)
    lite=(/966., 968., 970., 972., 974., 976., 978., 980., 982., 984., &
     & 986., 988., 990., 992., 994., 996., 998.,1000.,1005.,1010., &
     & 1015.,1020.,1025.,1030.,1035.,1040.,1045.,1050.,1055.,1060., &
     & 1065.,1070.,1075.,1080.,1085.,1090.,1095.,1100.,1105.,1110., &
     & 1115.,1120.,1125.,1130.,1135.,1140.,1145.,1150.,1155.,1160., &
     & 1165.,1170.,1175.,1180.,1185.,1190.,1195.,1200.,1205.,1210., &
     & 1215.,1220.,1225.,1230.,1235.,1240.,1245.,1250.,1255.,1260., &
     & 1265.,1270.,1275.,1280.,1285.,1290.,1295.,1300.,1305.,1310., &
     & 1315.,1320.,1325.,1330.,1335.,1340.,1345.,1350.,1355.,1360., &
     & 1365.,1370.,1375.,1380.,1385.,1390.,1395.,1400.,1405.,1410., &
     & 1415.,1420.,1425.,1430.,1435.,1440.,1445.,1450.,1455.,1458.81670, &
     & 1465.,1470.,1475.,1480.,1485.,1490.,1495.,1500.,1505.,1510., &
     & 1515.,1520.,1525.,1530.,1535.,1540.,1545.,1550.,1555.,1560., &
     & 1565.,1570.,1575.,1580.,1585.,1590.,1595.,1600.,1610.,1620., &
     & 1630.,1640.,1650.,1660.,1670.,1680.,1690.,1700.,1710.,1720., &
     & 1730.,1740.,1750.,1760.,1770.,1780.,1790.,1800.,1810.,1820., &
     & 1830.,1840.,1850.,1860.,1870.,1880.,1890.,1900.,1910.,1920., &
     & 1930.,1940.,1950.,1960.,1970.,1980.,1990.,2000.,2010.,2020., &
     & 2030.,2040.,2050.,2060.,2070.,2080.,2090.,2100.,2110.,2120./)

    litf=(/ 2130.,2140.,2150.,2160.,2170.,2180.,2190.,2200.,2210.,2220., &
     & 2230.,2240.,2250.,2260.,2270.,2279.40330,2290.,2300.,2310.,2320., &
     & 2330.,2340.,2350.,2360.,2370.,2380.,2390.,2400.,2410.,2420., &
     & 2430.,2440.,2450.,2460.,2470.,2480.,2490.,2500.,2510.,2520., &
     & 2530.,2540.,2550.,2560.,2570.,2580.,2590.,2600.,2610.,2620., &
     & 2630.,2640.,2650.,2660.,2670.,2680.,2690.,2700.,2710.,2720., &
     & 2730.,2740.,2750.,2760.,2770.,2780.,2790.,2800.,2810.,2820., &
     & 2830.,2840.,2850.,2860.,2870.,2880.,2890.,2900.,2910.,2920., &
     & 2930.,2940.,2950.,2960.,2970.,2980.,2990.,3000.,3010.,3020., &
     & 3030.,3040.,3050.,3060.,3070.,3080.,3090.,3100.,3110.,3120., &
     & 3130.,3140.,3150.,3160.,3170.,3180.,3190.,3200.,3220.,3240., &
     & 3260.,3282.34820,3300.,3320.,3340.,3360.,3380.,3400.,3420.,3440., &
     & 3460.,3480.,3500.,3520.,3540.,3560.,3580.,3600.,3620.,3640., &
     & 3660.,3680.,3700.,3720.,3740.,3760.,3780.,3800.,3820.,3840., &
     & 3860.,3880.,3900.,3920.,3940.,3960.,3980.,4000.,4020.,4040., &
     & 4060.,4080.,4100.,4120.,4140.,4160.,4180.,4200.,4220.,4240., &
     & 4260.,4280.,4300.,4320.,4340.,4360.,4380.,4400.,4420.,4440., &
     & 4460.,4480.,4500.,4520.,4540.,4560.,4580.,4600.,4620.,4640., &
     & 4660.,4680.,4700.,4720.,4740.,4760.,4780.,4800.,4820.,4840./)
    litg=(/4860.,4880.,4900.,4920.,4940.,4960.,4980.,5000.,5020.,5040., &
     & 5060.,5080.,5100.,5120.,5140.,5160.,5180.,5200.,5220.,5240., &
     & 5260.,5280.,5300.,5320.,5340.,5360.,5380.,5400.,5420.,5440., &
     & 5460.,5480.,5500.,5520.,5540.,5560.,5580.,5600.,5620.,5640., &
     & 5660.,5680.,5700.,5720.,5740.,5760.,5780.,5800.,5820.,5840., &
     & 5860.,5880.,5900.,5920.,5940.,5960.,5980.,6000.,6020.,6040., &
     & 6060.,6080.,6100.,6120.,6140.,6160.,6180.,6200.,6220.,6240., &
     & 6260.,6280.,6300.,6320.,6340.,6360.,6380.,6400.,6440.,6480., &
     & 6520.,6560.,6600.,6640.,6680.,6720.,6760.,6800.,6840.,6880., &
     & 6920.,6960.,7000.,7040.,7080.,7120.,7160.,7200.,7240.,7280., &
     & 7320.,7360.,7400.,7440.,7480.,7520.,7560.,7600.,7640.,7680., &
     & 7720.,7760.,7800.,7840.,7880.,7920.,7960.,8000.,8040.,8080., &
     & 8120.,8160.,8200.,8240.,8280.,8320.,8360.,8400.,8440.,8480., &
     & 8520.,8560.,8600.,8640.,8680.,8720.,8760.,8800.,8840.,8880., &
     & 8920.,8960.,9000.,9040.,9080.,9120.,9160.,9200.,9240.,9280., &
     & 9320.,9360.,9400.,9440.,9480.,9520.,9560.,9600.,9640.,9680., &
     & 9720.,9760.,9800.,9840.,9880.,9920.,9960.,10000./)

      wavebig(1:95)=biga
      wavebig(96:220)=bigb
      wavebig(221:nsizebig+1) = bigc
      wavelit(1:95)=lita
      wavelit(96:285)=litb
      wavelit(286:475)=litc
      wavelit(476:665)=litd
      wavelit(666:855)=lite
      wavelit(856:1045)=litf
      wavelit(1046:nsizelit+1)=litg

   else

      open(unit = 95, file = './INPUT/bin_grid_sizes.dat', form= 'formatted', &
     &     status = 'old', access = 'sequential' )
!
      read(95,*) nsizebig
      do i = 1, nsizebig
        read(95,*) wavebig(i)
      end do
      wavelit(1) = wavebig(1)
      wavelit(2) = wavebig(2)
      nsizebig = nsizebig -1
!
      close(unit = 95)

   end if

END SUBROUTINE def_binsize

END MODULE DFCOMM


text:
MODULE PRMTRS_STMAS_CLOUD

  ! LAPS PARAMETERS
  CHARACTER*9          :: LAPSAST            ! LAPS ASCII TIME
  INTEGER              :: LAPSI4T            ! LAPS I4TIME
  INTEGER              :: ICYCLE             ! LAPS CYCLE TIME

  REAL                 :: RMISSING           ! LAPS DEFAULT MISSING VALUE (YUANFU)
  REAL                 :: BADSFCDT           ! LAPS DEFAULT BAD SURFACE (YUANFU)
  REAL    ,ALLOCATABLE :: LATITUDE(:,:)      ! LAPS LATITUDE OF THE GRID (YUANFU)
  REAL    ,ALLOCATABLE :: LONGITUD(:,:)      ! LAPS LONGITUDE OF THE GRID (YUANFU)
  REAL    ,ALLOCATABLE :: TOPOGRPH(:,:)      ! LAPS LONGITUDE OF THE GRID (YUANFU)
  ! REAL    ,ALLOCATABLE :: LAPSRADW(:)        ! LAPS RADIAL WIND (YUANFU)
  REAL    ,ALLOCATABLE :: Z_FCSTGD(:)        ! VERTICAL LEVEL OF MODEL OR LAPS

  ! END OF LAPS PARAMETERS.

  ! GPS PROCESSING PARAMETERS:
  ! SPECIFIC GAS CONSTANT FOR WATER VAPOR:
  REAL    ,  PARAMETER :: RV_GAS=461
  ! COEFFICIENTS OF REFRACTIVITY AT MICROWAVE FREQUENCIES:
  REAL    ,  PARAMETER :: K2_RFC=22.1
  REAL    ,  PARAMETER :: K3_RFC=373900

  ! BUFR PARAMETERS:
  REAL    ,  PARAMETER :: BUFRMISS = 10.0E10 ! MISSING VALUE IN REAL FOR BUFR DATA
  REAL    ,  PARAMETER :: BADFLAG = -999.9   ! MISSING VALUE

  INTEGER ,  PARAMETER :: MAXDIMS=4          ! NUMBER OF TOTAL DIMENSION
  INTEGER ,  PARAMETER :: TMGOBS_CHANNEL=32  ! CHANNEL FOR TMG FILE
  INTEGER ,  PARAMETER :: PIGOBS_CHANNEL=33  ! CHANNEL FOR TMG FILE
  INTEGER              :: BGTHOBS            ! NOT USE BACKGROUND UNTIL NOBSGRID LESS THAN BGTHOBS
  INTEGER              :: GRDLEVL            ! GRID LEVEL
  INTEGER              :: NUMGRID(MAXDIMS)   ! NUMBER OF GRID FOR EVERY DIMENSION
  INTEGER              :: INIGRID(MAXDIMS)   ! INITIAL NUMBER OF GRID FOR EVERY DIMENSION
  INTEGER              :: MAXGRID(MAXDIMS)   ! MAX NUMBER OF GRID FOR EVERY DIMENSION
  INTEGER              :: NTMPGRD(MAXDIMS)   ! NUMBER OF OLD GRID FOR EVERY DIMENSION
  INTEGER              :: FCSTGRD(MAXDIMS)   ! NUMBER OF OLD GRID FOR EVERY DIMENSION
  INTEGER              :: NUMDIMS            ! NUMBER OF DIMENSION VALID
  INTEGER              :: NGPTOBS            ! NUMBER OF GRID POINT PER OBSERVATION
  INTEGER              :: NUMSTAT            ! NUMBER OF STATE
  INTEGER              :: NALLOBS            ! NUMBER OF TOTAL OBSERVATION VALID
  INTEGER ,ALLOCATABLE :: NOBSTAT(:)         ! NUMBER OF VARIETY OF OBSERVATION VALID
  INTEGER              :: NUMVARS            ! NUMBER OF TOTAL CONTROL VARIABLES VALID
  INTEGER              :: IFBKGND            ! WHETHER USE BACKGROUND OR NOT, 1 IS YES, 0 IS NOT
  INTEGER              :: IFBOUND            ! WHETHER USE BOUND OR NOT, 1 IS YES, 0 IS NOT
  INTEGER              :: IFPCDNT            ! WHETHER USE PURE PRESSURE COORDINATE, 1 IS YES, 0 IS NOT
  INTEGER              :: IF_TEST            ! WHETHER RUN FOR TEST CASE, 1 IS FOR THE TEST CASE, ADDED BY ZHONGJIE HE
  INTEGER              :: FNSTGRD            ! THE FINEST GRID LEVEL
  INTEGER              :: IFREPET            ! WHETHERE RUN THE MUTIGRID FRAME IN MONOTONOUS OR REPEATEDLY, 0 FOR MONOTONOUS, 1 FOR REPEATEDLY
  INTEGER              :: ITREPET            ! THE TIMES TO REPEAT

  REAL    ,ALLOCATABLE :: PENAL0X(:)         ! INITIAL X DIRECTION PENALTY COEFFICENT
  REAL    ,ALLOCATABLE :: PENAL0Y(:)         ! INITIAL Y DIRECTION PENALTY COEFFICENT
  REAL    ,ALLOCATABLE :: PENAL0Z(:)         ! INITIAL Z DIRECTION PENALTY COEFFICENT
  REAL    ,ALLOCATABLE :: PENAL0T(:)         ! INITIAL T DIRECTION PENALTY COEFFICENT
  REAL    ,ALLOCATABLE :: PENAL_X(:)         ! X DIRECTION PENALTY COEFFICENT
  REAL    ,ALLOCATABLE :: PENAL_Y(:)         ! Y DIRECTION PENALTY COEFFICENT
  REAL    ,ALLOCATABLE :: PENAL_Z(:)         ! Z DIRECTION PENALTY COEFFICENT
  REAL    ,ALLOCATABLE :: PENAL_T(:)         ! T DIRECTION PENALTY COEFFICENT
  REAL                 :: PNLT0PU            ! INITIAL GEOSTROPHIC BALANCE PENALTY COEFFICENT, FOR P AND U
  REAL                 :: PNLT0PV            ! INITIAL GEOSTROPHIC BALANCE PENALTY COEFFICENT, FOR P AND V
  REAL                 :: PNLT_PU            ! GEOSTROPHIC BALANCE PENALTY COEFFICENT, FOR P AND U
  REAL                 :: PNLT_PV            ! GEOSTROPHIC BALANCE PENALTY COEFFICENT, FOR P AND V
  REAL                 :: GRDSPAC(MAXDIMS)   ! GRID SPACING
  REAL                 :: ORIPSTN(MAXDIMS)   ! ORIGINAL POSITION OF THE STUDY DOMAIN
  REAL                 :: COSTFUN            ! COST FUNCTION
  REAL    ,ALLOCATABLE :: GRADINT(:,:,:,:,:) ! GRADIENT OF COST FUNCTION
  CHARACTER(LEN=100)   :: OBSFILE            ! OBSERVATION FILE NAME

  INTEGER              :: COSSTEP            ! THE ITERATE STEPS BEFORE MIDDLE GRID LEVEL
  INTEGER              :: MIDGRID            ! MIDDLE GRID LEVEL WHERE THE ITERATE STEP CHANGED
  INTEGER              :: FINSTEP            ! THE ITERATE STEPS AFTER MIDDLE GRID LEVEL

  INTEGER              :: U_CMPNNT           ! U COMPONENT INDEX
  INTEGER              :: V_CMPNNT           ! V COMPONENT INDEX
  INTEGER              :: W_CMPNNT           ! W COMPONENT INDEX
  INTEGER              :: PRESSURE           ! PRESSURE INDEX
  INTEGER              :: TEMPRTUR           ! TEMPERATURE INDEX
  INTEGER              :: HUMIDITY           ! SPECIFIC HUMIDITY INDEX
!added by shuyuan 20100721 for radar reflectivity
  INTEGER              :: ROUR_CMPNNT        ! DENSITY*RAIN WATER MIXING RATIO
  INTEGER              :: ROUS_CMPNNT        ! DENSITY*SNOW WATER MIXING RATIO
!
  INTEGER              :: XSL                ! X COORDINATE INDEX
  INTEGER              :: YSL                ! Y COORDINATE INDEX
  INTEGER              :: PSL                ! PRESSURE INDEX
  INTEGER              :: CSL                ! CORIOLIS FORCE INDEX
  INTEGER              :: DSL                ! DENSITY INDEX

  ! INTEGER ,ALLOCATABLE :: IDXRADWN(:,:)      ! GRID INDICES OF RADIAL WIND

  REAL                 :: XYTRANS            ! COEFFICIENT USED TO TRANSLATE THE X AND Y COORDINATE TO METERS
  REAL                 :: Z_TRANS            ! COEFFICIENT USED TO TRANSLATE THE Z COORDINATE TO METERS

  REAL    ,ALLOCATABLE :: WWW(:,:,:,:)       ! W COMPONENT
  REAL    ,ALLOCATABLE :: XXX(:,:)           ! DISTANCE OF HORIZONTAL DIMENSION 1
  REAL    ,ALLOCATABLE :: YYY(:,:)           ! DISTANCE OF HORIZONTAL DIMENSION 2
  REAL    ,ALLOCATABLE :: ZZZ(:,:,:,:)       ! HEIGHT OF GENERAL VERTICAL COORDINATE
  REAL    ,ALLOCATABLE :: PPP(:)             ! PRESSURE VERTICAL COORDINATE
  REAL    ,ALLOCATABLE :: COR(:,:)           ! CORIOLIS FREQUENCY
  REAL    ,ALLOCATABLE :: DEG(:,:)           ! DIRECTION ANGLE
  REAL    ,ALLOCATABLE :: DEN(:,:,:,:)       ! DENSITY
  REAL    ,ALLOCATABLE :: SCL(:)             ! SCALING THE PHYSICAL VARIABLE
  REAL    ,ALLOCATABLE :: SL0(:)             ! DEFAULT SCALING THE PHYSICAL VARIABLE
  REAL                 :: SCP(5)             ! SCALING XXX,YYY,ZZZ OR PPP,COR,DEN
  REAL                 :: ORIVTCL
  REAL    ,ALLOCATABLE :: XX0(:,:)           ! DISTANCE OF HORIZONTAL DIMENSION 1
  REAL    ,ALLOCATABLE :: YY0(:,:)           ! DISTANCE OF HORIZONTAL DIMENSION 2
  REAL    ,ALLOCATABLE :: ZZ0(:,:,:,:)       ! HEIGHT OF GENERAL VERTICAL COORDINATE
  REAL    ,ALLOCATABLE :: ZZB(:)       ! SIGMA LEVEL FOR VERTICAL COORDINATE
  REAL    ,ALLOCATABLE :: PP0(:)             ! PRESSURE VERTICAL COORDINATE
  REAL    ,ALLOCATABLE :: PPM(:)             ! MULTIGRID PRESSURE VERTICAL COORDINATE
  REAL    ,ALLOCATABLE :: CR0(:,:)           ! CORIOLIS FREQUENCY
  REAL    ,ALLOCATABLE :: DG0(:,:)           ! DIRECTION ANGLE
  REAL    ,ALLOCATABLE :: DN0(:,:,:,:)       ! DENSITY
! FOR OBSERVATION
  REAL                 :: OBSRADAR           ! THE OBSERVATION RATIO OF CONVENTION TO RADAR DATA
  REAL                 :: OBS_SFMR           ! THE OBSERVATION RATIO OF CONVENTION TO SFMR DATA
  REAL                 :: OBSREF           ! THE OBSERVATION RATIO OF REF
  REAL    ,ALLOCATABLE :: OBSPOSTN(:,:)      ! COORDINATE OF OBSERVATION POINT
  REAL    ,ALLOCATABLE :: OBSCOEFF(:,:)      ! INTERPOLATION COEFFICENT
  INTEGER ,ALLOCATABLE :: OBSIDXPC(:,:)      ! INDEX OF EVERY NODE PER CELL
  INTEGER ,ALLOCATABLE :: OBSSTATE(:)        ! NUMBER OF OBSERVATION FOR EVERY STATE
  REAL    ,ALLOCATABLE :: OBSVALUE(:)        ! OBSERVATION VALUE
  REAL    ,ALLOCATABLE :: OBSERROR(:)        ! OBSERVATION ERROR
  REAL    ,ALLOCATABLE :: OBSEINF1(:)        ! OBSERVATION EXTRA INFORMATION 1
  REAL    ,ALLOCATABLE :: OBSEINF2(:)        ! OBSERVATION EXTRA INFORMATION 2
  REAL    ,ALLOCATABLE :: OBSEINF3(:)        ! OBSERVATION EXTRA INFORMATION 3

  ! GPSMET WETDELAY DATA:
  INTEGER, PARAMETER :: MAX_GPS = 5000		! TEMPORARILY HARD CODED  !!changefrom 2000 to 5000 20100525 liu
  INTEGER :: NUM_GPS,I4T_GPS			! I4T_GPS: I4TIME OF Jan 1, 1970
  REAL    :: GPS_XYT(3,MAX_GPS),GPS_ELV(MAX_GPS),GPS_TIM(MAX_GPS), &
             GPS_ERR(MAX_GPS),GPS_TPW(MAX_GPS),GPS_WET(MAX_GPS)

! FOR GRID POINT VARIABLE
  REAL    ,ALLOCATABLE :: GRDBKGD0(:,:,:,:,:)! GRID BACKGROUND VALUE
  REAL    ,ALLOCATABLE :: GRDBKGND(:,:,:,:,:)! GRID BACKGROUND VALUE
  REAL    ,ALLOCATABLE :: GRDANALS(:,:,:,:,:)! GRID ANALYSIS VALUE
  REAL    ,ALLOCATABLE :: TMPANALS(:,:,:,:,:)! GRID ANALYSIS VALUE

  REAL    ,ALLOCATABLE :: WWWOUT(:,:,:,:)    ! W COMPONENT FOR OUTPUT
  REAL    ,ALLOCATABLE :: DIFFTOUT(:,:,:,:,:)! GRID ANALYSIS VALUE FOR OUTPUR
  REAL    ,ALLOCATABLE :: Z_MAXGID(:)        ! VERTICAL LEVEL OF FINEST ANALYSIS GRID

  REAL                 :: PNLT0HY            ! INITIAL HYDROSTATIC CONDITION PENALTY COEFFICENT. ADDIED BY ZHONGJIE HE
  REAL                 :: PNLT_HY            ! HYDROSTATIC CONDITION PENALTY COEFFICENT. ADDIED BY ZHONGJIE HE
  REAL                 :: TAUL_HY            ! REDUCTION COEFFICENT OF HYDROSTATIC CONDITION PENALTY TERM, ADDIED BY ZHONGJIE HE
  REAL                 :: ENDHYLV            ! THE GRID LEVEL AFTER WHICH THE HYDROSTATIC CONDITION PENALTY TERM IS OMMITTED. ADDED BY ZHONGJIE HE
  REAL                 :: ENDGSLV            ! THE GRID LEVEL AFTER WHICH THE GEOSTROPHIC BALANCE PENALTY TERM IS OMMITTED. ADDED BY ZHONGJIE HE
  REAL                 :: LIMIT_3            ! THE LIMITATION OF HEIGHT TO DECIDE IN WHICH RANGE THE OBSERVATION IS AVIABLE. ADDED BY ZHONGJIE HE
  REAL                 :: LIMIT_4            ! THE LIMITATION OF TIME TO DECIDE IN WHICH RANGE THE OBSERVATION IS AVIABLE. ADDED BY ZHONGJIE HE

  ! ANALYSIS DOMAIN BY YUANFU XIE:
  INTEGER              :: ITIME2(2)          ! ANALYSIS TIME WINDOW

!RUIXIA LIU ADD 20120321
  REAL    ,ALLOCATABLE :: AMSUA(:,:,:,:)
  REAL    ,ALLOCATABLE :: T_SFC_0(:,:,:)
  REAL    ,ALLOCATABLE :: AMSUAGRD(:,:,:,:)  ! SATELLITE MULTIGRID POSITION
  REAL    ,ALLOCATABLE :: T_SFC_GRD(:,:,:)
  REAL    ,ALLOCATABLE :: P_lyr(:)
  REAL    ,ALLOCATABLE :: ozn_lyr(:)
  REAL    ,ALLOCATABLE :: P_lvl(:)
  REAL    ,ALLOCATABLE :: T_SFC(:,:)
  REAL    ,ALLOCATABLE :: PPP1(:)             ! PRESSURE VERTICAL COORDINATE
  REAL    ,ALLOCATABLE :: rad(:,:,:,:)
  REAL    ,ALLOCATABLE :: ja(:,:,:,:,:)
  REAL    ,ALLOCATABLE :: T_lyr(:,:,:)
  REAL    ,ALLOCATABLE :: Q_lyr(:,:,:)
  REAL  :: CMAG
  INTEGER  :: INUM2
  CHARACTER(LEN=100)   :: SAT_DA_USE
  CHARACTER(LEN=200)   :: Co_Path
  CHARACTER(LEN=100)   :: Sat_Data_Path(2)   ! ADD AMSU-B BY YUANFU FEB 2013
  CHARACTER(LEN=100)   :: Moist_Unit(2)   ! ADD FOR AMSU-A/B BY YUANFU FEB 2013
  INTEGER  :: NCHBEG,NCHEND
  INTEGER  :: NAMSUB ! YUANFU ADDED FOR COUNTING NUMBER AMSU-B OBS
  INTEGER  ,ALLOCATABLE  :: AMSUID(:,:)
  INTEGER  :: AMSUA_CHANNELS_USED(15),AMSUB_CHANNELS_USED(5) ! YUANFU ADDED FOR CONTROLING AMSU DATA
  REAL  ,ALLOCATABLE :: T_SFC_BKG0(:,:,:)
  REAL    ,ALLOCATABLE :: T_lyr1(:,:,:)
  REAL    ,ALLOCATABLE :: Q_lyr1(:,:,:)
  REAL    ,ALLOCATABLE ::  T_SFC1(:,:)
  REAL    ,ALLOCATABLE :: PPP2(:)
  REAL OZN_lyr1(41)
  REAL  ,  ALLOCATABLE :: rad1(:,:,:,:)
  REAL    ,ALLOCATABLE :: ja1(:,:,:,:,:)
  REAL  ,  ALLOCATABLE :: RAD_AMSUB(:,:,:,:)
  REAL    ,ALLOCATABLE :: JAC_AMSUB(:,:,:,:,:)
  REAL    ,ALLOCATABLE :: P_lvl1(:)
  REAL    ,ALLOCATABLE :: LND(:,:)
  REAL    ,ALLOCATABLE :: LND0(:,:)
  REAL    ,ALLOCATABLE :: LNDGRD(:,:)
  REAL    ,ALLOCATABLE :: LANDTYPE(:)
  REAL    ,ALLOCATABLE :: LNDTY0(:)
!END RUIXIA LIU ADD 20120321

  ! YUANFU XIE ADDS DATA STRUCTURES FOR SATELLITE DATA:
  TYPE SATELLITE_OBS
    INTEGER :: NUMPRO,NUMCHS  ! NUMBER OF PROFILES AND CHANNELS
    INTEGER, ALLOCATABLE :: GRDIDX(:,:)    ! HORIZONTAL GRID LOCATIONS
    INTEGER, ALLOCATABLE :: LNDCVR(:)      ! LAND COVERAGE: LAND:1 WATER:0
    INTEGER, ALLOCATABLE :: LNDTYP(:)      ! LAND TYPE: CRTM LAND OR WATER TYPE
    REAL, ALLOCATABLE    :: SATVAL(:,:)    ! OBSERVED SATELLITE VALUES: CHANNELS X PROFILES
    REAL, ALLOCATABLE    :: ANGLES(:,:)    ! SCAN & ZENITH ANGLES: 4 X PROFILES
                                           ! (LOCAL & SOLAR ANGLES)
  END TYPE SATELLITE_OBS

  TYPE(SATELLITE_OBS) :: AMSUB(3)          ! YUANFU ADDED AMSUB FEB 2013, MAX 3 TIME FRAME NOW
  TYPE(SATELLITE_OBS) :: AMSUBM(3)         ! YUANFU ADDED AMSUB FOR A MULTIGRID LEVEL FEB 2013, MAX 3 TIME FRAME NOW

CONTAINS

!=================================================================================
!>
!! THIS ROUTINE INITIALIES A SATELLITE OBSERVATION DATA STRUCTURE.
!! 1. ALLOCATE MEMORY;
!! 2. ASSIGN VALUES.
!!
!! \history \n
!! CREATION: MAR. 2013 \n
!! AUTHOR  : YUANFU XIE AT NOAA/ESRL/GSD
!
!=================================================================================
SUBROUTINE SATOBS_INIT(SAT_OBS,NUM_PRO,NUM_CHS,GRD_IDX,LND_CVR,LND_TYP,SAT_VAL,SAT_ANG)

  IMPLICIT NONE

  INTEGER, INTENT(IN) :: NUM_PRO,NUM_CHS,GRD_IDX(2,NUM_PRO)
  INTEGER, INTENT(IN) :: LND_CVR(NUM_PRO),LND_TYP(NUM_PRO)
  REAL,    INTENT(IN) :: SAT_VAL(NUM_CHS,NUM_PRO),SAT_ANG(4,NUM_PRO)

  TYPE(SATELLITE_OBS) :: SAT_OBS

  ! LOCAL VARIABLES:
  INTEGER :: IP,ISTATUS

  ! ALLOCATE:
  ALLOCATE(SAT_OBS%GRDIDX(2,NUM_PRO),SAT_OBS%LNDCVR(NUM_PRO),SAT_OBS%LNDTYP(NUM_PRO), &
           SAT_OBS%SATVAL(NUM_CHS,NUM_PRO),SAT_OBS%ANGLES(4,NUM_PRO), STAT=ISTATUS)

  SAT_OBS%NUMPRO = NUM_PRO
  SAT_OBS%NUMCHS = NUM_CHS

  DO IP=1,NUM_PRO
    SAT_OBS%GRDIDX(:,IP) = GRD_IDX(:,IP)
    SAT_OBS%LNDCVR(IP) = LND_CVR(IP)
    SAT_OBS%LNDTYP(IP) = LND_TYP(IP)
    SAT_OBS%SATVAL(:,IP) = SAT_VAL(:,IP)
    SAT_OBS%ANGLES(:,IP) = SAT_ANG(:,IP)
  ENDDO

END SUBROUTINE SATOBS_INIT

!=================================================================================
!>
!! THIS ROUTINE PROJECTS A SATELLITE OBSERVATION DATA STRUCTURE ON TO ONE OVER A
!! FINER RESOLUTION.
!!
!! \HISTORY \N
!! CREATION: MAR. 2013
!!
!! \AUTHOR: YUANFU XIE AT NOAA/ESRL/GSD
!
!=================================================================================
SUBROUTINE SATOBS_PROJ(SATOBS1,SATOBS2,LEVEL,IFLAG)

  IMPLICIT NONE

  INTEGER, INTENT(IN)  :: LEVEL
  INTEGER, INTENT(OUT) :: IFLAG
  TYPE(SATELLITE_OBS), INTENT(IN)  :: SATOBS1
  TYPE(SATELLITE_OBS), INTENT(OUT) :: SATOBS2

    ! LOCAL VARIABLES:
  INTEGER :: IP,ISTATUS,IBASE
  REAL    :: SQRT

  IF (LEVEL .LE. 0) THEN
    IFLAG = 0
    PRINT*,'NO PROJECTION IS NEEDED'
    RETURN
  ENDIF
  IBASE = 2**(LEVEL-1)

  SATOBS2%NUMCHS = SATOBS1%NUMCHS
  SATOBS2%NUMPRO = SATOBS1%NUMPRO ! WORST CASE SCENARIO: PROJECTION HAS ALL OBS

  ! ALLOCATE:
  ALLOCATE(SATOBS2%GRDIDX(2,SATOBS2%NUMPRO),SATOBS2%LNDCVR(SATOBS2%NUMPRO), &
           SATOBS2%LNDTYP(SATOBS2%NUMPRO),SATOBS2%ANGLES(4,SATOBS2%NUMPRO), &
           SATOBS2%SATVAL(SATOBS2%NUMCHS,SATOBS2%NUMPRO), STAT=ISTATUS)

  SATOBS2%NUMPRO = 0
  DO IP=1,SATOBS1%NUMPRO
    IF (MOD(SATOBS1%GRDIDX(1,IP)-1,IBASE) .EQ. 0 .AND. &
        MOD(SATOBS1%GRDIDX(2,IP)-1,IBASE) .EQ. 0) THEN
      SATOBS2%NUMPRO = SATOBS2%NUMPRO + 1

      SATOBS2%GRDIDX(:,SATOBS2%NUMPRO) = (SATOBS1%GRDIDX(:,IP)-1)/(2**(LEVEL-1))+1
      SATOBS2%LNDCVR(SATOBS2%NUMPRO) = SATOBS1%LNDCVR(IP)
      SATOBS2%LNDTYP(SATOBS2%NUMPRO) = SATOBS1%LNDTYP(IP)
      SATOBS2%SATVAL(:,SATOBS2%NUMPRO) = SATOBS1%SATVAL(:,IP)
      SATOBS2%ANGLES(:,SATOBS2%NUMPRO) = SATOBS1%ANGLES(:,IP)
    ENDIF
  ENDDO

  IFLAG = 1

END SUBROUTINE SATOBS_PROJ

!=================================================================================
!>
!! THIS	FOLLOWING ROUTINE IS TO DEALLOCATE A SATELLITE_OBS DATA TYPE VARIABLE.
!!
!! \HISTORY \N
!! CREATION: MAR. 2013 \N
!! \AUTHOR: YUANFU XIE AT NOAA/ESRL
!
!=================================================================================

SUBROUTINE SATOBS_DEAD(SATOBS)

  TYPE(SATELLITE_OBS) :: SATOBS

  ! LOCAL VARIABLES:
  INTEGER :: ISTATUS

  SATOBS%NUMPRO = 0
  SATOBS%NUMCHS = 0

  DEALLOCATE(SATOBS%GRDIDX,SATOBS%LNDCVR,SATOBS%LNDTYP,SATOBS%SATVAL,SATOBS%ANGLES, &
             STAT=ISTATUS)

END SUBROUTINE SATOBS_DEAD

END MODULE PRMTRS_STMAS_CLOUD

max_stars_repo_name: maxinye/laps-mirror

max_stars_repo_licenses: Intel Unlicense OLDAP-2.2.1 NetCDF

max_issues_repo_name: longwosion/laps-mirror

max_issues_repo_licenses: Intel NetCDF OLDAP-2.2.1 Unlicense

max_forks_repo_name: longwosion/laps-mirror

max_forks_repo_licenses: Intel NetCDF OLDAP-2.2.1 Unlicense


text:
      SUBROUTINE GRF_CURVE
     &          (xarray,yarray,nrpts,xmin,ymin,xscale,yscale,
     &           line,dots_per_inch)

c***********************************************************************

c  Function:
c     Draws a curve representing the input x,y pairs

c  Parameters passed:
c     xarray        [r] array of X values; user units
c     yarray        [r] array of Y values; user units
c     nrpts         [i] number of points to be plotted
c     xmin          [r] origin of the X values; user units
c     ymin          [r] origin of the Y values; user units
c     xscale        [r] scale along the X axis; user units per inch
c     yscale        [r] scale along the Y axis; user units per inch
c     line          [i] line type to be used while drawing the curve
c     dots_per_inch [i] number of dost per inch; used to define the
c                       length of segments of all line types

c     NOTE: line=1: solid
c                2: long dash
c                3: medium dash
c                4: short dash
c                5: dotted
c                6: short + long dash
c                7: short + short + long dash

c  Parameters returned:

c  Common blocks referenced:

c  Functions and subroutines referenced:
c     mod
c     sqrt

c     grf_draw
c     grf_move

c  Entry:

c  References:

c  Change History:

c*******************!***************************************************

      include      'graphics.cmn'

      integer       dots_per_inch

      dimension     xarray(1),yarray(1)

      logical       up(8),move
      integer       jo(7),jc(7)

      data          up/.false.,.false.,.false.,.true.,
     &                 .false.,.true. ,.false.,.true./,
     &              jo/1,1,2,3,4,1,1/,
     &              jc/1,4,3,2,1,6,8/

c     If insufficient number of points, don't do anything
      if (nrpts .le. 1) RETURN

      if (graphics_device(1:3) .eq. 'sys') then

c        Set up screen output

c        Set up parameters for line type
         resolution=1./dots_per_inch

         kk=MOD(line,7)
         if (kk .le. 0) kk=kk+7
         j_origin=jo(kk)
         j_cycle =jc(kk)
         move=up(j_origin)
         j=0

c        Two different integrals are formed:
c        The first one is the total length of the curve being drawn <rho>;
c        the second is the length of the curve with respect to a segment of
c        the particular line type being used <rho2>.
         rho1=0.
         rho2=resolution
         dr=0.

c        Select the first point
         px1=(xarray(1)-xmin)/xscale
         py1=(yarray(1)-ymin)/yscale

c        Go to first position with pen up
         call GRF_MOVE (px1,py1)

         if (kk .eq. 5) then

c           The curve is dotted; generate a very short dash along the line
c           segment between the first and second points
            px2=(xarray(2)-xmin)/xscale
            py2=(yarray(2)-ymin)/yscale
            delx=px2-px1
            dely=py2-py1
            rho=SQRT(delx**2+dely**2)
            if (rho .gt. 0.)
     &         call GRF_DRAW (px1+.01*delx/rho,py1+.01*dely/rho)
         end if

c        Start curve
         do i=2,nrpts

c           Select next point
            px2=(xarray(i)-xmin)/xscale
            py2=(yarray(i)-ymin)/yscale

            if (kk .eq. 1) then

c              Solid line:
c              No need to integrate the curve
               call GRF_DRAW (px2,py2)

            else

c              Multi-segmented line:
c              Integrate the curve to get even-length line segments

c              First, calculate the length of the segment
               delx=px2-px1
               dely=py2-py1
               rho=SQRT(delx**2+dely**2)

c              Increment integration variable
               rho1=rho1+rho
               if (rho2 .le. rho1) then

c                 The data point is beyond the end of the line segment.
                  delx=delx*resolution/rho
                  dely=dely*resolution/rho

                  if (dr .gt. 0.) then

c                    Complete current line segment
                     dx=delx*dr/resolution
                     dy=dely*dr/resolution
                     px1=px1+dx
                     py1=py1+dy

                     if (kk .eq. 5) then
c                       Dotted line
                        call GRF_MOVE (px1,py1)
                        call GRF_DRAW (px1+.1*delx,py1+.1*dely)
                     else
     &               if (move) then
                        call GRF_MOVE (px1,py1)
                     else
                        call GRF_DRAW (px1,py1)
                     end if

c                    Increment interval counter
                     j=j+1
                     move=up(j_origin+MOD(j,j_cycle))
                     rho2=rho2+resolution
                  end if

                  do while (rho2 .le. rho1)
                     px1=px1+delx
                     py1=py1+dely

                     if (kk .eq. 5) then
c                       Dotted line
                        call GRF_MOVE (px1,py1)
                        call GRF_DRAW (px1+.1*delx,py1+.1*dely)
                     else
     &               if (move) then
                        call GRF_MOVE (px1,py1)
                     else
                        call GRF_DRAW (px1,py1)
                     end if

c                    Increment interval counter
                     j=j+1
                     move=up(j_origin+MOD(j,j_cycle))
                     rho2=rho2+resolution
                  end do
               end if

c              Complete segment
               if (move) then
                  call GRF_MOVE (px2,py2)
               else
                  call GRF_DRAW (px2,py2)
               end if
            end if
            dr=rho2-rho1
            px1=px2
            py1=py2
         end do
      else

c        Get the pattern length as a percentage of the diagonal
c        distance of the plot page
         pl=400./dots_per_inch/
     &      SQRT((graphics_extent(1)/graphics_units(1))**2+
     &           (graphics_extent(2)/graphics_units(2))**2)

c        Set the line type
         kk=MOD(line,7)
         if (kk .le. 0) kk=kk+7
         if (kk .eq. 1) then
            write(graphics_LogicalUnit,'('' lt;'')')
         else
            write(graphics_LogicalUnit,'('' lt'',i1,f6.2'';'')')
     &            kk-1,pl
         end if

c        Select the first point
         px1=(xarray(1)-xmin)/xscale
         py1=(yarray(1)-ymin)/yscale

c        Go to first position with pen up
         call GRF_MOVE (px1,py1)

c        Start curve
         do i=2,nrpts

c           Select next point
            px2=(xarray(i)-xmin)/xscale
            py2=(yarray(i)-ymin)/yscale

            call GRF_DRAW (px2,py2)
         end do

c        Reset line type to solid line
         write(graphics_LogicalUnit,'('' lt;'')')
      end if

      RETURN
      END      ! GRF_CURVE

ext: for

max_stars_repo_name: spinorkit/LWPC

max_issues_repo_name: spinorkit/LWPC

max_forks_repo_name: spinorkit/LWPC


text:
module modis19

	use CKD_common

	implicit none
	private

	integer, parameter :: ni = 4

	public :: ch19, init_ch19

contains

	subroutine init_ch19(modis_channel)

		type(ckd_type), intent(inout) :: modis_channel
		integer :: all_weights

		all_weights = ni

		modis_channel%num_weights = all_weights
		allocate(modis_channel%weights(all_weights))
		allocate(modis_channel%taus(all_weights, max_layers))

	end subroutine init_ch19

      subroutine ch19(u0, p0, t0, ux, fac, modis_channel)
!
!      This routine was developed for MODIS Channel 19
!                        10430--10960 cm^{-1}
!
!-----------------------------------------------------------------
!     Channel 19   0.912-0.959 microns
!
!     revised by David Kratz on 7/26/2007
!     also accounts for the continuum (routine modis19x.f90 does not)
!
!-----------------------------------------------------------------
!  INPUTS:
!     u0 --> water vapor amount [g/cm2/km]
!     p0 --> pressure in atmospheres
!     t0 --> temp in K
!     ux --> ozone [g/cm2/km]
!     dz --> layer thickness in km
!-----------------------------------------------------------------

		real, dimension(:), intent(in) :: u0, p0, t0, ux, fac
		type(ckd_type), intent(inout) :: modis_channel

		integer, parameter :: nlev = max_levels
		integer, parameter :: nlay = max_layers
		integer, parameter :: mlv = nlay

		real :: p(nlay), t(nlay),u(nlay), &
				uo2(nlay),uo3(nlay)

		real :: tau(ni,nlay), f(ni), twinx(ni), twin

		integer :: m, i

      do  m=1,use_layers
        t(m)=(t0(m)+t0(m+1))/2.0
        p(m)=(p0(m)+p0(m+1))/2.0
        u(m)=fac(m)*(u0(m)+u0(m+1))/2.0
      enddo

      call ck_19(u,f,p,t,tau)

      do  m=1,use_layers
        call window_19(u(m),twin,p(m),t(m),fac(m))
        twinx(1)=0.1634974*twin*1.643519
        twinx(2)=0.3509185*twin*1.643519
        twinx(3)=0.3403536*twin*1.643519
        twinx(4)=0.1452305*twin*1.643519

         do i=1,ni ! 4 k's for h2o
          modis_channel%taus(i,m)=tau(i,m)+twinx(i)
 			if (m == 1) modis_channel%weights(i) = f(i)
         end do
      enddo

      end subroutine ch19

! *********************************************************************
      subroutine ck_19(u,f,p,t,tau)
		real, dimension(:), intent(in) :: p, t, u

		real, dimension(:), intent(inout) :: f
		real, dimension(:,:), intent(inout) :: tau

		real :: k(ni)
		real :: coefk(ni,3,num_pressures)
		integer :: i, jp, jt

      f(1)=0.537422
      f(2)=0.347788
      f(3)=0.101705
      f(4)=0.013085
      k(1)=0.040109846
      do i=2,ni
        k(i)=11.0*k(i-1)
      end do
      data ( (coefk(1,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.00125,0.00185,0.00273,0.00405,0.00601,0.00896, &
     & 0.01341,0.02018,0.03045,0.04619,0.06981,0.10605, &
     & 0.16153,0.24574,0.37054,0.55524,0.82447,1.20557, &
     &  1.70015, &
     &  5.125E-06, 7.000E-06, 9.250E-06, 1.262E-05, 1.737E-05, 2.362E-05, &
     &  3.250E-05, 4.562E-05, 6.550E-05, 9.500E-05, 1.390E-04, 2.065E-04, &
     &  3.081E-04, 4.695E-04, 7.158E-04, 1.017E-03, 1.432E-03, 2.051E-03, &
     &  2.800E-03, &
     &  3.125E-09,-6.250E-09,-1.250E-08,-2.813E-08,-4.063E-08,-6.563E-08, &
     & -1.063E-07,-1.844E-07,-2.750E-07,-4.500E-07,-7.125E-07,-1.019E-06, &
     & -1.491E-06,-2.012E-06,-2.431E-06,-2.937E-06,-4.234E-06,-5.688E-06, &
     & -8.038E-06/
      data ( (coefk(2,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.02466,0.02586,0.02771,0.03036,0.03419,0.03973, &
     & 0.04756,0.05880,0.07466,0.09742,0.12914,0.17471, &
     & 0.23979,0.33231,0.45988,0.63405,0.86191,1.14948, &
     &  1.47436, &
     &  3.060E-04, 3.040E-04, 3.026E-04, 3.030E-04, 3.066E-04, 3.120E-04, &
     &  3.224E-04, 3.389E-04, 3.619E-04, 3.895E-04, 4.204E-04, 4.584E-04, &
     &  5.066E-04, 5.650E-04, 6.431E-04, 6.987E-04, 7.815E-04, 9.155E-04, &
     &  1.079E-03, &
     &  1.156E-06, 1.131E-06, 1.066E-06, 1.019E-06, 9.781E-07, 9.000E-07, &
     &  8.656E-07, 7.469E-07, 7.156E-07, 5.938E-07, 3.656E-07, 2.344E-07, &
     &  9.688E-08,-1.875E-07,-4.219E-07,-1.025E-06,-1.862E-06,-3.950E-06, &
     & -4.584E-06/
      data ( (coefk(3,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.29496,0.29539,0.29641,0.29760,0.29960,0.30325, &
     & 0.30870,0.31819,0.33261,0.35556,0.38639,0.43257, &
     & 0.49728,0.58131,0.68510,0.80666,0.93435,1.06129, &
     &  1.16180, &
     &  1.732E-03, 1.730E-03, 1.721E-03, 1.712E-03, 1.702E-03, 1.675E-03, &
     &  1.641E-03, 1.588E-03, 1.513E-03, 1.405E-03, 1.256E-03, 1.092E-03, &
     &  9.172E-04, 7.450E-04, 5.989E-04, 4.624E-04, 3.570E-04, 1.972E-04, &
     &  4.012E-05, &
     & -1.287E-06,-1.269E-06,-1.356E-06,-1.316E-06,-1.147E-06,-1.206E-06, &
     & -9.031E-07,-9.906E-07,-6.344E-07,-7.406E-07,-7.344E-07,-5.875E-07, &
     & -1.125E-06,-1.019E-06,-1.459E-06,-2.522E-06,-2.181E-06,-1.275E-06, &
     & -1.416E-06/
      data ( (coefk(4,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 1.78024,1.77949,1.78023,1.77836,1.77539,1.77265, &
     & 1.76533,1.75779,1.74214,1.72111,1.67675,1.62309, &
     & 1.55539,1.47017,1.36110,1.22980,1.07615,0.91625, &
     &  0.75218, &
     & -1.576E-03,-1.573E-03,-1.594E-03,-1.586E-03,-1.551E-03,-1.556E-03, &
     & -1.528E-03,-1.487E-03,-1.429E-03,-1.370E-03,-1.275E-03,-1.171E-03, &
     & -1.064E-03,-9.365E-04,-8.485E-04,-7.939E-04,-7.222E-04,-6.522E-04, &
     & -5.941E-04, &
     &  1.056E-06, 1.047E-06, 4.562E-07, 4.312E-07, 9.875E-07, 3.687E-07, &
     &  1.506E-06, 2.250E-07, 1.319E-06, 7.688E-07, 1.501E-07, 7.750E-07, &
     &  6.250E-07, 6.250E-07, 1.694E-06, 2.278E-06, 2.963E-06, 2.788E-06, &
     &  1.797E-06/

		call apply_ckd(ni, k, u, f, p, t, tau, coefk, .false.)

      end subroutine ck_19
! *********************************************************************
      subroutine window_19 (amnt,twin,patm,temp,fac)
!     Parameterization of CKD_2.4 continuum over 10430 to 10960 cm-1 band
!     Fit inaccurate for amnt < 0.01 g/cm*2;
!          however optical depths < 0.001.
!     INPUT:
!     amnt = h2O layer amount (g/cm**2)
!     patm= pressure (atm)
!     temp = temperature (k)
!     OUTPUT:
!     twin = parameterized CKD_2.1 optical depth over band
!

	  real, intent(in) :: amnt, patm, temp, fac
	  real, intent(inout) :: twin

      integer, parameter :: ncoef=7

      real ::  aa(ncoef), ph2o, tau_log

      data aa/-2.119,1.002,-5.192e-03,0.9057,4.698,-1.894e-03,5.871e-02/

      ph2o = amnt*(8.314e+07*temp)/(fac*1.0e+05*18.01534*1.01325e+06)

          tau_log   = aa(1)              + &
                      aa(2) * log(amnt)  + &
                      aa(3) * temp       + &
                      aa(4) * patm       + &
                      aa(5) * ph2o       + &
                      aa(6) * amnt       + &
                      aa(7) * log(ph2o)

        twin = exp ( tau_log )

        end subroutine window_19

end module modis19

max_stars_repo_name: GEOS-ESM/AeroApps

max_issues_repo_name: GEOS-ESM/AeroApps

max_forks_repo_name: GEOS-ESM/AeroApps