text:
C    %W% %G%
      real function tzlead (t2, t1, dt, x0, x1, y0, ndivs)
      implicit none
      real t2, t1, dt, x0, x1, y0
      integer ndivs

c     -  solves a first order lead block:
c              sT2
c            -------
c            1 + sT1
c
c     -  Uses sub-time step logic if ndivs > 1
c
c     -   subscript 0 means value at start of interval
c             "     1  "      "   "   end  "     "
c     -   x is input var
c         y is output var
c
      real aa,bb,v0,v1,w0,w1
      real rdivs,r0wt,r1wt,va,vb,wa,wb,dts
      integer la

c     -     begin     begin     begin     begin     begin     begin

      if (t1 .eq. 0.0) then
        tzlead = x1
        return
      endif
      if (dt .eq. 0.0) then
        tzlead = y0 + ((x1 - x0) * t2 / t1)
        return
      endif

c     -  No time step subdivisions

      if (ndivs .le. 1) then
        aa = dt / 2.0 / t1
        bb = t2 * (x0 + x1) / t1
        v0 = (x0 * t2 / t1) - y0
        v1 = ( v0 - aa * v0 + aa * bb) / (aa + 1.0)
        tzlead = (x1 * t2 / t1) - v1
        return
      endif
c     -  Using time step subdivisions.  Find v1 = F (w0,w1,v0)
      rdivs = ndivs
      r0wt = rdivs
      r1wt = 0.0
      dts = dt / rdivs
      aa = dts / 2.0 / t1
      w0 = x0 * t2 / t1
      w1 = x1 * t2 / t1
      wa = w0
      va = w0 - y0
      do la = 1,ndivs
        r0wt = r0wt - 1.0
        r1wt = r1wt + 1.0
        wb = ((r0wt * w0) + (r1wt * w1)) / rdivs
        bb = wa + wb
        vb = ( va - aa * va + aa * bb) / (aa + 1.0)
        wa = wb
        va = vb
      enddo
      tzlead = wb - vb
c
      return
      end

max_stars_repo_name: mbheinen/bpa-ipf-tsp

max_issues_repo_name: cuihantao/bpa-ipf-tsp

max_forks_repo_name: mbheinen/bpa-ipf-tsp


text:
module probdata_module

  ! rad shock parameters
  use amrex_fort_module, only : rt => amrex_real
  real(rt)        , save ::  rho_0, T_0, rhoe_0

  real(rt)        , save :: xmin,xmax,ymin,ymax,zmin,zmax

end module probdata_module

max_stars_repo_name: yingtchen/Castro

max_issues_repo_name: yingtchen/Castro

max_forks_repo_name: yingtchen/Castro


text:
      program f77iterate_c
C
C    This example program illustrates how to use the CFITSIO iterator function.
C
C    This program creates a 2D histogram of the X and Y columns of an event
C    list.  The 'main' routine just creates the empty new image, then executes
C    the 'writehisto' work function by calling the CFITSIO iterator function.
C
C    'writehisto' opens the FITS event list that contains the X and Y columns.
C    It then calls a second work function, calchisto, (by recursively calling
C    the CFITSIO iterator function) which actually computes the 2D histogram.

C     external work function to be passed to the iterator
      external writehisto

      integer ncols
      parameter (ncols=1)
      integer units(ncols), colnum(ncols), datatype(ncols)
      integer iotype(ncols), offset, n_per_loop, status
      character*70 colname(ncols)

      integer naxes(2), ounit, blocksize
      character*80 fname
      logical exists

C     include f77.inc -------------------------------------
C     Codes for FITS extension types
      integer IMAGE_HDU, ASCII_TBL, BINARY_TBL
      parameter (
     &     IMAGE_HDU  = 0,
     &     ASCII_TBL  = 1,
     &     BINARY_TBL = 2  )

C     Codes for FITS table data types

      integer TBIT,TBYTE,TLOGICAL,TSTRING,TSHORT,TINT
      integer TFLOAT,TDOUBLE,TCOMPLEX,TDBLCOMPLEX
      parameter (
     &     TBIT        =   1,
     &     TBYTE       =  11,
     &     TLOGICAL    =  14,
     &     TSTRING     =  16,
     &     TSHORT      =  21,
     &     TINT        =  31,
     &     TFLOAT      =  42,
     &     TDOUBLE     =  82,
     &     TCOMPLEX    =  83,
     &     TDBLCOMPLEX = 163  )

C     Codes for iterator column types

      integer InputCol, InputOutputCol, OutputCol
      parameter (
     &     InputCol       = 0,
     &     InputOutputCol = 1,
     &     OutputCol      = 2  )
C     End of f77.inc -------------------------------------

C**********************************************************************
C     Need to make these variables available to the 2 work functions
      integer xsize,ysize,xbinsize,ybinsize
      common /histcomm/ xsize,ysize,xbinsize,ybinsize
C**********************************************************************

      status = 0

      xsize = 480
      ysize = 480
      xbinsize = 32
      ybinsize = 32

      fname = 'histoimg.fit'
      ounit = 15

C     delete previous version of the file if it exists
      inquire(file=fname,exist=exists)
      if( exists ) then
         open(ounit,file=fname,status='old')
         close(ounit,status='delete')
      endif
 99   blocksize = 2880

C     create new output image
      call ftinit(ounit,fname,blocksize,status)

      naxes(1) = xsize
      naxes(2) = ysize

C     create primary HDU
      call ftiimg(ounit,32,2,naxes,status)

      units(1) = ounit

C     Define column as TINT and Output
      datatype(1) = TINT
      iotype(1) = OutputCol

C     force whole array to be passed at one time
      n_per_loop = -1
      offset = 0

C     execute the function to create and write the 2D histogram
      print *,'Calling writehisto iterator work function... ',status

      call ftiter( ncols, units, colnum, colname, datatype, iotype,
     &      offset, n_per_loop, writehisto, 0, status )

      call ftclos(ounit, status)

C     print out error messages if problem
      if (status.ne.0) then
         call ftrprt('STDERR', status)
      else
        print *,'Program completed successfully.'
      endif

      stop
      end

C--------------------------------------------------------------------------
C
C   Sample iterator function.
C
C   Iterator work function that writes out the 2D histogram.
C   The histogram values are calculated by another work function, calchisto.
C
C--------------------------------------------------------------------------
      subroutine writehisto(totaln, offset, firstn, nvalues, narrays,
     &     units_out, colnum_out, datatype_out, iotype_out, repeat,
     &     status, userData, histogram )

      integer totaln,offset,firstn,nvalues,narrays,status
      integer units_out(narrays),colnum_out(narrays)
      integer datatype_out(narrays),iotype_out(narrays)
      integer repeat(narrays)
      integer histogram(*), userData

      external calchisto
      integer ncols
      parameter (ncols=2)
      integer units(ncols), colnum(ncols), datatype(ncols)
      integer iotype(ncols), rowoffset, rows_per_loop
      character*70 colname(ncols)

      integer iunit, blocksize
      character*80 fname

C     include f77.inc -------------------------------------
C     Codes for FITS extension types
      integer IMAGE_HDU, ASCII_TBL, BINARY_TBL
      parameter (
     &     IMAGE_HDU  = 0,
     &     ASCII_TBL  = 1,
     &     BINARY_TBL = 2  )

C     Codes for FITS table data types

      integer TBIT,TBYTE,TLOGICAL,TSTRING,TSHORT,TINT
      integer TFLOAT,TDOUBLE,TCOMPLEX,TDBLCOMPLEX
      parameter (
     &     TBIT        =   1,
     &     TBYTE       =  11,
     &     TLOGICAL    =  14,
     &     TSTRING     =  16,
     &     TSHORT      =  21,
     &     TINT        =  31,
     &     TFLOAT      =  42,
     &     TDOUBLE     =  82,
     &     TCOMPLEX    =  83,
     &     TDBLCOMPLEX = 163  )

C     Codes for iterator column types

      integer InputCol, InputOutputCol, OutputCol
      parameter (
     &     InputCol       = 0,
     &     InputOutputCol = 1,
     &     OutputCol      = 2  )
C     End of f77.inc -------------------------------------

C**********************************************************************
C     Need to make these variables available to the 2 work functions
      integer xsize,ysize,xbinsize,ybinsize
      common /histcomm/ xsize,ysize,xbinsize,ybinsize
C**********************************************************************

      if (status .ne. 0) return

C     name of FITS table
      fname = 'iter_c.fit'
      iunit = 16

C     do sanity checking of input values
      if (totaln .ne. nvalues) then
C     whole image must be passed at one time
         status = -1
         return
      endif

      if (narrays .ne. 1) then
C     number of images is incorrect
         status = -2
         return
      endif

      if (datatype_out(1) .ne. TINT) then
C     input array has wrong data type
         status = -3
         return
      endif

C     open the file and move to the table containing the X and Y columns
      call ftopen(iunit,fname,0,blocksize,status)
      call ftmnhd(iunit, BINARY_TBL, 'EVENTS', 0, status)
      if (status) return

C     both the columns are in the same FITS file
      units(1) = iunit
      units(2) = iunit

C     desired datatype for each column: TINT
      datatype(1) = TINT
      datatype(2) = TINT

C     names of the columns
      colname(1) = 'X'
      colname(2) = 'Y'

C     leave column numbers undefined
      colnum(1) = 0
      colnum(2) = 0

C     define whether columns are input, input/output, or output only
C     Both input
      iotype(1) = InputCol
      iotype(1) = InputCol

C     take default number of rows per iteration
      rows_per_loop = 0
      rowoffset = 0

C     calculate the histogram
      print *,'Calling calchisto iterator work function... ', status

      call ftiter( ncols, units, colnum, colname, datatype, iotype,
     &      rowoffset, rows_per_loop, calchisto, histogram, status )

      call ftclos(iunit,status)
      return
      end

C--------------------------------------------------------------------------
C
C   Iterator work function that calculates values for the 2D histogram.
C
C--------------------------------------------------------------------------
      subroutine calchisto(totalrows, offset, firstrow, nrows, ncols,
     &     units, colnum, datatype, iotype, repeat, status,
     &     histogram, xcol, ycol )

      integer totalrows,offset,firstrow,nrows,ncols,status
      integer units(ncols),colnum(ncols),datatype(ncols)
      integer iotype(ncols),repeat(ncols)
      integer histogram(*),xcol(*),ycol(*)
C     include f77.inc -------------------------------------
C     Codes for FITS extension types
      integer IMAGE_HDU, ASCII_TBL, BINARY_TBL
      parameter (
     &     IMAGE_HDU  = 0,
     &     ASCII_TBL  = 1,
     &     BINARY_TBL = 2  )

C     Codes for FITS table data types

      integer TBIT,TBYTE,TLOGICAL,TSTRING,TSHORT,TINT
      integer TFLOAT,TDOUBLE,TCOMPLEX,TDBLCOMPLEX
      parameter (
     &     TBIT        =   1,
     &     TBYTE       =  11,
     &     TLOGICAL    =  14,
     &     TSTRING     =  16,
     &     TSHORT      =  21,
     &     TINT        =  31,
     &     TFLOAT      =  42,
     &     TDOUBLE     =  82,
     &     TCOMPLEX    =  83,
     &     TDBLCOMPLEX = 163  )

C     Codes for iterator column types

      integer InputCol, InputOutputCol, OutputCol
      parameter (
     &     InputCol       = 0,
     &     InputOutputCol = 1,
     &     OutputCol      = 2  )
C     End of f77.inc -------------------------------------

      integer ii, ihisto, xbin, ybin

C**********************************************************************
C     Need to make these variables available to the 2 work functions
      integer xsize,ysize,xbinsize,ybinsize
      common /histcomm/ xsize,ysize,xbinsize,ybinsize
C**********************************************************************

      if (status .ne. 0) return

C    --------------------------------------------------------
C      Initialization procedures: execute on the first call
C    --------------------------------------------------------
      if (firstrow .eq. 1) then
C     do sanity checking of input values

         if (ncols .ne. 2) then
C     number of arrays is incorrect
            status = -4
            return
         endif

         if (datatype(1).ne.TINT .or. datatype(2).ne.TINT) then
C     wrong datatypes
            status = -5
            return
         endif

C     initialize the histogram image pixels = 0, including null value
         do 10 ii = 1, xsize * ysize + 1
            histogram(ii) = 0
 10     continue

      endif

C     ------------------------------------------------------------------
C       Main loop: increment the 2D histogram at position of each event
C     ------------------------------------------------------------------

      do 20 ii=2,nrows+1
        xbin = xcol(ii) / xbinsize
        ybin = ycol(ii) / ybinsize

        ihisto = ( ybin * xsize ) + xbin + 2
        histogram(ihisto) = histogram(ihisto) + 1
 20   continue

      return
      end

ext: f

max_stars_repo_name: tyronedong/OpenGR

max_issues_repo_name: tyronedong/OpenGR

max_forks_repo_name: tyronedong/OpenGR


text:
c
c*****************************
c
        subroutine rdcndwt(w,psiprime)

        parameter (u0=2.8)

        if (w.gt. 4) then
        w=0.0
        psiprime=0.0
        return
        end if

        t=(-exp(u0*(w-u0)))
        u=w
        w=exp(t)
        t=u0*t*u

        if((t+1.).lt.(1.e-20)) then
           psiprime=0.
           return
        end if

        psiprime=(w*(1.+t))

        return
        end

max_stars_repo_name: simpeg-research/iris-mt-scratch

max_issues_repo_name: simpeg-research/iris-mt-scratch

max_forks_repo_name: simpeg-research/iris-mt-scratch


text:
module linear_solvers
  !! Module for computing solutions to linear systems of the form:
  !! $$\mathbf{A}\mathbf{x}=\mathbf{b}$$
  implicit none

  contains

  function solve_U(U,b) result(x)
    !! Computes the solution of an upper triangular system:
    !! $$ \mathbf{U}\mathbf{x} = \mathbf{b}$$
      real, intent(in) :: U(:,:)
    !! Upper triangular matrix
      real, intent(in) :: b(:)
    !! Right hand side vector
      real :: x(size(U,2))
    !! Solution vector **x**

    integer :: i, N

    N = size(U,2)

    x(N) = b(N)/U(N,N)

    do i=N-1,1,-1
      x(i) = ( b(i) - dot_product(U(i,i+1:N),x(i+1:N)))/U(i,i)
    end do

  end function

  function solve_L(L,b) result(x)
  !! Computes the solution of a lower triangular system:
  !! $$ \mathbf{L}\mathbf{x} = \mathbf{b}$$
    real, intent(in) :: L(:,:)
  !! Lower triangular matrix
    real, intent(in) ::  b(:)
  !! Right hand side vector
    real :: x(size(L,2))
  !! Solution vector **x**

    integer :: i, N

    N = size(L,2)

    x(1) = b(1)/L(1,1)

    do i=2,N
      x(i) = ( b(i) - dot_product(L(i,1:i-1),x(1:i-1)))/L(i,i)
    end do

  end function

  subroutine make_U(M)
    !! Takes the extended matrix M = (A|b)
    !! and makes it upper triangular.
    real, intent(inout) :: M(:,:)
      !! Extended matrix (A|b)
    integer ::  N, i, j, imax
    real :: Q

    real :: temp(size(M,2))

    N = size(M,1)
    do j = 1, N-1
      ! row swapping for stability
      imax = MAXLOC(abs(M(j:N,j)),1) + (j - 1)
      temp = M(j,:)
      M(j,:) = M(imax,:)
      M(imax,:) = temp
      ! triangularization of column j
      do i = j+1, N
          Q = M(i,j)/M(j,j)
          M(i,:) = M(i,:) - Q*M(j,:)
      end do
    end do

  end subroutine

  function gauss_solve(A,b) result(x)
    !! Computes the solution of the system:
    !!  $$\mathbf{A}\mathbf{x}=\mathbf{b}$$
    !! where **A** is a non-singular matrix.
      real, intent(in) :: A(:,:)
    !! Non-singular matrix
      real, intent(in) :: b(:)
    !! Right hand side vector
      real :: x(size(A,2))
    !! Solution vector **x**

    real :: M(size(A,1),size(A,2)+1)
    integer :: N

    N = size(A,2)
    M(:,1:N) = A
    M(:,N+1) = b

    call make_U(M)

    x = solve_U(M(:,1:N),M(:,N+1))
  end function

    subroutine factor_LU(A,L,U)
      !! Computes de LU factorization of a non-singular matrix **A**
      !! $$\mathbf{A}=\mathbf{L}\mathbf{U}$$
        real, intent(in) :: A(:,:)
      !! Non-singular matrix
        real, intent(out) :: L(size(A,1),size(A,2))
      !! Lower triangular matrix
        real, intent(out) :: U(size(A,1),size(A,2))
      !! Upper triangular matrix

      integer :: N, i, j, k

      N = size(A,2)

      L = 0
      U = 0

      do k = 1, N

        do j = k, N
          U(k,j) = A(k,j) - dot_product(L(k,1:k-1),U(1:k-1,j))
        end do

        L(k,k) = 1
        do i = k+1, N
          L(i,k) = ( A(i,k) - dot_product(L(i,1:k-1),U(1:k-1,k))) / U(k,k)
        end do

      end do

    end subroutine

    function solve_LU(L,U,b) result(x)
      !! Computes the solution to the system:
      !! $$\mathbf{L}\mathbf{U}\mathbf{x}=\mathbf{b}$$
      !! where **L** and **U** are the matrices obtained from
      !! a LU factorization.
        real, intent(in) :: L(:,:)
      !! Lower triangular matrix
      real, intent(in) :: U(:,:)
      !! Upper triangular matrix
        real, intent(in) :: b(:)
      !! Right hand side vector
        real :: x(size(b))
      !! Solution vector **x**

      x = solve_U(U,solve_L(L,b))

    end function

    function solve(A,b) result(x)
      !! Computes the solution to the system:
      !! $$\mathbf{A}\mathbf{x}=\mathbf{b}$$
      !! using LU factorization.
        real, intent(in) :: A(:,:)
      !! Non-singular matrix
        real, intent(in) :: b(:)
      !! Right hand side vector
        real :: x(size(b))
      !! Solution vector **x**

      real :: L(size(A,1),size(A,2)), U(size(A,1), size(A,2)), Ap(size(A,1),size(A,2))

      call factor_LU(A,L,U)

      x = solve_LU(L,U,b)

    end function

    function Jacobi(A,b,x_0, maxIter, err_x) result(x)
      !! Computes the solution to the system
      !! $$\mathbf{A}\mathbf{x}=\mathbf{b}$$
      !! using the Jacobi iterative method.
        real, intent(in) :: A(:,:)
      !! Non-singular matrix
        real, intent(in) :: b(:)
      !! Right hand side vector
        real, intent(in) :: x_0(:)
      !! Initial approximation for the iteration
        integer, intent(in) :: maxIter
      !! Maximum number of iterations
        real, intent(in) :: err_x
      !! Precision for the stop criterion
        real :: x(size(A,2))
      !! Solution vector **x**

      integer :: i, N
      real :: x_old(size(A,2)), L(size(A,1),size(A,2)),U(size(A,1),size(A,2)), D_inv(size(A,2))

      N = size(A,2)

      L = 0
      U = 0
      do i = 1, N
        L(i,1:i-1)=A(i,1:i-1)
        D_inv(i) = 1./A(i,i)
        U(i,i+1:N) = A(i,i+1:N)
      end do

      x_old = x_0
      do i = 1, maxIter
        x = -D_inv*MATMUL(L+U,x_old)+D_inv*b
        if ( norm2(x-x_old)/norm2(x) < err_x ) exit
        x_old = x
      end do
      if( i>maxIter ) print*, 'the maximum number of iteration was reached &
                              & without obtaining convergence'
    end function

    function GaussSeidel(A,b,x_0, maxIter, err_x) result(x)
      !! Computes the solution to the system
      !! $$\mathbf{A}\mathbf{x}=\mathbf{b}$$
      !! using the Gauss-Seidel iterative method.
        real, intent(in) :: A(:,:)
      !! Non-singular matrix
        real, intent(in) :: b(:)
      !! Right hand side vector
        real, intent(in) :: x_0(:)
      !! Initial approximation for the iteration
        integer, intent(in) :: maxIter
      !! Maximum number of iterations
        real, intent(in) :: err_x
      !! Precision for the stop criterion
        real :: x(size(A,2))
      !! Solution vector **x**

      integer :: i, N
      real :: x_old(size(A,2)), LD(size(A,1),size(A,2)), U(size(A,1),size(A,2))

      N = size(A,2)

      LD = 0
      U = 0
      do i = 1, N
        LD(i,  1:i) = A(i,  1:i)
        U (i,i+1:N) = A(i,i+1:N)
      end do

      x_old = x_0
      do i = 1, maxIter
        x = solve_L(LD,-MATMUL(U,x_old)+b)
        if ( norm2(x-x_old)/norm2(x) < err_x ) exit
        x_old = x
      end do
      if( i>maxIter ) print*, 'the maximum number of iteration was reached &
                              & without obtaining convergence'
    end function

end module

max_stars_repo_name: MPenaR/NumericalMethods

max_issues_repo_name: MPenaR/NumericalMethods

max_forks_repo_name: MPenaR/NumericalMethods


! testcov.f90 --
!     Program to instrument Fortran code (free form) so that the
!     test coverage can be determined and reported.
!
!     Note:
!     The type of coverage this program provides is "branch coverage",
!     a fairly basic form of test coverage.
!
program testcov
    use testcov_reg

    implicit none

    character(len=200) :: line
    character(len=200) :: inname
    character(len=200) :: outname

    integer            :: ierr
    logical            :: report
    logical            :: exists

    !
    ! Statement types
    !
    integer, parameter :: is_comment            = 1
    integer, parameter :: is_declaration        = 2
    integer, parameter :: is_other              = 3
    integer, parameter :: is_endsubprogram      = 4
    integer, parameter :: is_endprogram         = 5
    integer, parameter :: is_subprogram         = 6
    integer, parameter :: is_program            = 7
    integer, parameter :: if_statement          = 8
    integer, parameter :: else_statement        = 9
    integer, parameter :: endif_statement       = 10
    integer, parameter :: do_statement          = 11
    integer, parameter :: cycle_statement       = 12
    integer, parameter :: exit_statement        = 13
    integer, parameter :: enddo_statement       = 14
    integer, parameter :: select_statement      = 15
    integer, parameter :: case_statement        = 16
    integer, parameter :: casedefault_statement = 17
    integer, parameter :: endselect_statement   = 18
    integer, parameter :: goto_statement        = 19
    integer, parameter :: continue_statement    = 20

    open( 10, file = 'testcov.inp', status = 'old' )

    report = .false.

    do
        read( 10, '(a)', iostat = ierr ) line
        if ( ierr /= 0 ) then
            exit
        endif

        line = adjustl(line)

        if ( line(1:1) == '!' ) then
            cycle
        endif

        if ( line == 'MODE=report' ) then
            call testcov_load
            report = .true.
            cycle
        endif
        if ( line == 'MODE=instrument' ) then
            report = .false.
            cycle
        endif

        !
        ! Any other line is considered the name of a source file
        !
        inname = line
        inquire( file = inname, exist = exists )
        if ( .not. exists ) then
            write(*,*) 'File "'// trim(inname) //'" does not exist - skipping'
            cycle
        endif

        write(*,*) 'Processing "'// trim(inname) //'"'

        call make_outname( inname, outname, report )

        open( 20, file = inname )
        open( 21, file = outname )

        call process_file( report )

        close( 20 )
        close( 21 )
    enddo

    write(*,*) 'Done'

contains

! make_outname --
!     Determine the name of the output file
!
! Arguments:
!     inname         Name of the input file
!     outname        Name of the output file
!     report         Report coverage or not (instrument instead)
!
! Result:
!     If the input is "mysource.f90" then the result is either
!     "_mysource.f90" (for instrumenting) or "mysource.lst" (for
!     reporting)
!
subroutine make_outname( inname, outname, report )
    character( len=*), intent(in)  :: inname
    character( len=*), intent(out) :: outname
    logical, intent(in)            :: report

    integer                        :: k
    character(len=2)               :: slash = '/\\' ! Some compilers see \ as an escape character ...

    if ( report ) then
        k = index( inname, '.', .true. )
        if ( k == 0 ) then
            k = len_trim(line) ! Highly unlikely ..
        endif

        outname = inname(1:k) // 'lst'
    else
        k = index( inname, slash(1:1), .true. )
        if ( k == 0 ) then
            k = index( inname, slash(2:2), .true. )
        endif
        outname = inname(1:k-1) // '_' // inname(k+1:)
    endif
end subroutine make_outname

! process_file --
!     Read the source code and produce an instrumented file or a
!     coverage report
!
! Arguments:
!     report         Report coverage or not (instrument instead)
!
subroutine process_file( report )
    logical, intent(in)            :: report

    integer                        :: ierr
    integer                        :: k
    character(len=2)               :: slash = '/\\' ! Some compilers see \ as an escape character ...
    character(len=200), dimension(40) :: line

    logical                        :: start_sub
    logical                        :: start_program
    integer                        :: lineno
    integer                        :: codetype
    integer                        :: nolines
    integer                        :: i

    lineno        = 0
    start_sub     = .false.
    start_program = .true.    ! The case of a missing program statement!
    do
        !
        ! TODO: continuation lines!
        !
        call read_line( line, nolines, ierr )
        if ( ierr /= 0 ) then
            exit
        endif

        lineno = lineno + nolines

        call classify_line( line(1), codetype )
       !write( 21, '(i5,a)' ) codetype, trim(line(1))

        select case( codetype )
            case ( is_comment, is_declaration )
                if ( report ) then
                    write( 21, '(10x,a)' ) (trim(line(i)) ,i=1,nolines)
                else
                    write( 21, '(a)' ) (trim(line(i)) ,i=1,nolines)
                endif
                cycle

            case ( is_subprogram, is_program )
                if ( codetype == is_program ) then
                    start_program = .true.
                endif
                start_sub = .true.
                if ( report ) then
                    write( 21, '(10x,a)' ) (trim(line(i)) ,i=1,nolines)
                else
                    write( 21, '(a)' ) (trim(line(i)) ,i=1,nolines)
                endif
                cycle

            case ( is_endprogram, is_endsubprogram )
                if ( start_program ) then
                    start_program = .false.
                    call handle_start( line, lineno, nolines, codetype )
                else
                    if ( report ) then
                        write( 21, '(10x,a)' ) (trim(line(i)) ,i=1,nolines)
                    else
                        write( 21, '(a)' ) (trim(line(i)) ,i=1,nolines)
                    endif
                endif
                cycle

            case ( is_other )
                if ( start_sub ) then
                    start_sub = .false.
                    call handle_start( line, lineno, nolines, codetype )
                else
                    if ( report ) then
                        write( 21, '(10x,a)' ) (trim(line(i)) ,i=1,nolines)
                    else
                        write( 21, '(a)' ) (trim(line(i)) ,i=1,nolines)
                    endif
                endif
                cycle

            case ( if_statement, else_statement, endif_statement )
                if ( start_sub ) then
                    start_sub = .false.
                    call handle_start( line, lineno, nolines, codetype )
                endif
                call handle_if( line, lineno, nolines, codetype )

            case ( do_statement, cycle_statement, exit_statement, enddo_statement )
                if ( start_sub ) then
                    start_sub = .false.
                    call handle_start( line, lineno, nolines, codetype )
                endif
                call handle_do( line, lineno, nolines, codetype )

            case ( select_statement, case_statement, casedefault_statement, endselect_statement )
                if ( start_sub ) then
                    start_sub = .false.
                    call handle_start( line, lineno, nolines, codetype )
                endif
                call handle_case( line, lineno, nolines, codetype )

            case ( goto_statement, continue_statement )
                if ( start_sub ) then
                    start_sub = .false.
                    call handle_start( line, lineno, nolines, codetype )
                endif
!               call handle_goto( line, lineno, codetype )

            case default
                write(*,*) 'Programming error: unknown statement type!'
                write(*,*) 'Line: ', trim(line(1))
                write( 21, '(a)' ) (trim(line(i)) ,i=1,nolines)

        end select
    enddo
end subroutine process_file

! tolower --
!     Convert a string to lower case
!
! Arguments:
!     string         String to be converted
!
! Result:
!     String in lowercase
!
function tolower( string )
    character(len=*), intent(in) :: string
    character(len=len(string))   :: tolower

    character(len=26), parameter :: upper = 'ABCDEFGHIJKLMNOPQRSTUVWXYZ'
    character(len=26), parameter :: lower = 'abcdefghijklmnopqrstuvwxyz'

    integer                    :: i
    integer                    :: k
    integer                    :: length

    length  = len(string)
    tolower = string
    do i = 1,length
        k = index( upper, string(i:i) )
        if ( k > 0 ) then
           tolower(i:i) = lower(k:k)
        endif
    enddo
end function tolower

! is_digit --
!     Is a character a digit or not
!
! Arguments:
!     char           Character to be examined
!
! Result:
!     True if digit, false otherwise
!
logical function is_digit( char )
    character(len=1), intent(in) :: char

    is_digit = index( '01234567890', char ) > 0

end function is_digit

! is_letter --
!     Is a character a (lower/uppercase) letter or not
!
! Arguments:
!     char           Character to be examined
!
! Result:
!     True if letter, false otherwise
!
logical function is_letter( char )
    character(len=1), intent(in) :: char

    character(len=52), parameter :: letters = 'ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz'

    is_letter = index( letters, char ) > 0

end function is_letter

! adjust_line --
!     Convert to lower case and remove superfluous spaces and tabs
!
! Arguments:
!     line           Line to be adjusted
!
subroutine adjust_line( line )
    character(len=*), intent(inout)  :: line

    integer                          :: i
    integer                          :: j

    line = adjustl( tolower(line) )

    do i =1,len_trim(line)
        if ( line(i:i) == char(9) ) then
            line(i:i) = ' '
        endif
    enddo

    j = 0
    do i = 1,len_trim(line)
        if ( line(i:i) /= ' ' ) then
            j = j + 1
            line(j:j) = line(i:i)
        elseif ( line(i-1:i-1) /= ' ' ) then
            j = j + 1
            line(j:j) = line(i:i)
        endif
    enddo

    line(j+1:) = ' '

end subroutine adjust_line

! classify_line --
!     Determine the class of the line of code
!
! Arguments:
!     line           Line of code to examine
!     codetype       Class to which it belongs
!
subroutine classify_line( line, codetype )
    character(len=*), intent(in)    :: line
    integer, intent(out)            :: codetype

    character(len=200)              :: adjusted_line
    character(len=200)              :: next_piece

    integer                         :: i
    integer                         :: idxkey
    integer                         :: k
    logical                         :: found
    character(len=1)                :: type_next
    character(len=20)               :: keyword

    integer, dimension(10)          :: subprogram_type = &
    (/ is_endsubprogram, is_endprogram, &
       is_endprogram, is_endsubprogram, &
       is_endsubprogram, is_endsubprogram, &
       is_endsubprogram, &
       is_program, is_subprogram, is_subprogram /)
    character(len=20), dimension(10) :: subprogram_keyword = &
    (/ 'contains', 'end program', 'endprogram', 'end subroutine', &
       'endsubroutine', 'end function', 'endfunction', &
       'program', 'subroutine', 'function' /)

    ! Note: "end" is treated here, as we need to check the entire line

    integer, dimension(15)           :: control_type = &
    (/ do_statement, enddo_statement, enddo_statement, &
       if_statement, else_statement, if_statement, if_statement, &
       endif_statement, endif_statement, &
       select_statement, endselect_statement, endselect_statement, &
       casedefault_statement, case_statement, is_endsubprogram /)
    character(len=20), dimension(15) :: control_keyword = &
    (/ '$do', '%enddo', '%end do', '#if', '@else', '#elseif', '#else if', &
       '%endif', '%end if', '#select case', '%endselect', '%end select',  &
       '@case default', '#case', '@end' /)

    character(len=20), dimension(71) :: declare_keyword = &
    (/ '#real::',     '#real ::',     '#real,',     '#real ,',     '#real(',     '#real (',     '%real',     &
       '#integer::',  '#integer ::',  '#integer,',  '#integer ,',  '#integer(',  '#integer (',  '%integer',  &
       '#logical::',  '#logical ::',  '#logical,',  '#logical ,',  '#logical(',  '#logical (',  '%logical',  &
       '#complex::',  '#complex ::',  '#complex,',  '#complex ,',  '#complex(',  '#complex (',  '%complex',  &
       '#character::','#character ::','#character,','#character ,','#character(','#character (','%character',&
       '#double precision::', '#double precision,', '#double precision ,', '%double precision',              &
       '#type(', '#type (', '%type', '%end type', '%endtype',                                                &
       '#implicit none', '%implicit',                                                                        &
       '%use', '#data (', '%data',                                                                           &
       '#parameter::', '#parameter ::', '#parameter,', '#parameter ,', '#parameter(', '#parameter (',        &
       '#save::',      '#save ::',      '%save',                                                             &
       '#private::',   '#private ::',   '%private',                                                          &
       '#public::',    '#public ::',    '%public',                                                           &
       '%interface',   '%end interface','%endinterface',                                                     &
       '#external::',  '#external ::',  '%external',                                                         &
       '%module procedure'                                                                                &
    /)
    ! target, pointer, allocatable, module ...
    ! Oh, and of course: subroutine, function inside interface block

    !
    ! Massage the line we read: lower case characters and no multiple spaces
    !
    codetype = is_other
    adjusted_line = line
    call adjust_line( adjusted_line )
    found = .false.

    !
    ! Comments
    !
    if ( adjusted_line(1:1) == '!' .or. adjusted_line == ' ' ) then
        codetype = is_comment
        return
    endif

    !
    ! Program, subroutine or function
    !
    do i = 1,size(subprogram_keyword)
        k = index( adjusted_line, trim(subprogram_keyword(i)) // ' ' )
        if ( k > 0 ) then
            k = k + len_trim(subprogram_keyword(i))
            next_piece = adjustl( adjusted_line(k:) )
            if ( is_letter(next_piece(1:1)) .or. next_piece(1:1) == ' ') then
                codetype = subprogram_type(i)
                found = .true.
                exit
            endif
        endif
    enddo

    if ( found ) return

    !
    ! Do, if, case ...
    !
    do i = 1,size(control_keyword)
        type_next = control_keyword(i)(1:1)
        keyword   = control_keyword(i)(2:)
        idxkey    = i

        k = index( adjusted_line, trim(keyword) )
        if ( k > 0 ) then
            select case ( type_next )
                case ( '@' )
                    if ( adjusted_line == keyword ) then
                        found = .true.
                        exit
                    endif
                case ( '$' )
                    k = k + len_trim(keyword) + 1
                    if ( adjusted_line(k-1:k-1) == ' ' .and. &
                         ( is_letter(adjusted_line(k:k)) .or. &
                           is_digit(adjusted_line(k:k))  .or. &
                           adjusted_line == keyword ) )       then
                        found = .true.
                        exit
                    endif
                case ( '#' )
                    k = k + len_trim(keyword)
                    if ( adjusted_line(k:k)     == '(' .or. &
                         adjusted_line(k+1:k+1) == '('      ) then
                        found = .true.
                        exit
                    endif
                case ( '%' )
                    k = k + len_trim(keyword)
                    if ( adjusted_line(k:k) == ' ' ) then
                        found = .true.
                        exit
                    endif
            end select
        endif
    enddo

    if ( found ) then
        codetype = control_type(idxkey)
        return
    endif

    !
    ! Declarations ...
    !
    do i = 1,size(declare_keyword)
        type_next = declare_keyword(i)(1:1)
        keyword   = declare_keyword(i)(2:)
        idxkey    = i

        k = index( adjusted_line, trim(keyword) )
        if ( k == 1 ) then
            select case ( type_next )
                case ( '#' )
                    k = k + len_trim(keyword) + 1
                    if ( is_letter(adjusted_line(k:k)) .or. &
                         adjusted_line(k:k) == '_' ) then
                        found = .true.
                        exit
                    endif
                case ( '%' )
                    k = k + len_trim(keyword) + 2
                    if ( is_letter(adjusted_line(k:k)) .or. &
                         adjusted_line(k:k) == '_' ) then
                        found = .true.
                        exit
                    endif
            end select
        endif
    enddo

    if ( found ) then
        codetype = is_declaration
        return
    endif

end subroutine classify_line

! read_line --
!     Read the line and possible continuation lines
!
! Arguments:
!     line              Line of code
!     nolines           Number of lines
!     ierr              Error code (if any)
!
subroutine read_line( line, nolines, ierr )
    character(len=*), dimension(:), intent(inout) :: line
    integer, intent(out)                          :: nolines
    integer, intent(out)                          :: ierr

    integer                                       :: i
    integer                                       :: length

    nolines = size(line)
    do i = 1,size(line)
        read( 20, '(a)', iostat = ierr ) line(i)
        length = len_trim(line(i))
        if ( line(i)(length:length) /= '&' ) then
            nolines = i
            exit
        endif
    enddo

end subroutine read_line

! write_line --
!     Write the intrumentation code (low-level)
!
! Arguments:
!     line              Line of code
!     lineno            Line number
!     nolines           Number of lines
!     inname            Name of source file
!     code              What registration code to use
!
subroutine write_line( line, lineno, nolines, inname, code )
    character(len=*), dimension(:), intent(in) :: line
    integer, intent(in)                        :: lineno
    integer, intent(in)                        :: nolines
    character(len=*), intent(in)               :: inname
    integer, intent(in)                        :: code

    integer                                    :: i

    write( 21, '(a)' ) (trim(line(i)) ,i=1,nolines)
    write( 21, '(3a,i0,a,i0,a)' ) &
        '        call testcov_register__("', trim(inname), '",', &
        lineno, ',', code, ')'
end subroutine write_line

! handle_start --
!     Output for the start of the subprogram
!
! Arguments:
!     line           Line of code to instrument/report
!     lineno         Line number
!     nolines        Number of lines
!     codetype       Class to which it belongs
!
subroutine handle_start( line, lineno, nolines, codetype )
    character(len=*), dimension(:), intent(in)    :: line
    integer, intent(in)                           :: lineno
    integer, intent(in)                           :: nolines
    integer, intent(in)                           :: codetype

    integer                                       :: count
    integer                                       :: i

    if ( report ) then
        if ( codetype == is_endprogram .or. codetype == is_endsubprogram ) then
            write( 21, '(10x,a)' ) trim(line(1))
        else
            call testcov_get_count( inname, lineno, 1, count )
            write( 21, '(i8,a,a)' ) count, ': ', trim(line(1))
            if ( count == 0 ) then
                 write( 21, '(a8,a)'  ) '==> ', ': (subroutine/function never called)'
             endif
             if ( nolines > 1 ) then
                 write( 21, '(10x,a)'  ) ( trim(line(i)), i = 2,nolines )
             endif
         endif
    else
        !
        ! If the routine starts with an IF or DO or SELECT, then
        ! the number of entries will be counted anyway
        !
        if ( codetype == is_other ) then
            call write_line( (/ ' ' /), lineno, 1, inname, 1 )
        endif
        if ( codetype == is_endprogram .or. codetype == is_endsubprogram ) then
            write( 21, '(a)' ) '      call testcov_register__("",0,0)'
            write( 21, '(a)' ) ( trim(line(i)), i = 1,nolines )
        endif
    endif
end subroutine handle_start

! handle_do --
!     Output for the control statements of a do-loop
!
! Arguments:
!     line           Line of code to instrument/report
!     lineno         Line number
!     nolines        Number of lines
!     codetype       Class to which it belongs
!
! Note:
!     FORTRAN 77 style do-loops not treated correctedly!
!
subroutine handle_do( line, lineno, nolines, codetype )
    character(len=*), dimension(:), intent(in)    :: line
    integer, intent(in)                           :: lineno
    integer, intent(in)                           :: nolines
    integer, intent(in)                           :: codetype

    integer                                       :: i
    integer, dimension(100), save                 :: do_start = 0
    integer, save                                 :: do_level = 0
    integer                                       :: count
    integer                                       :: skipped

    if ( report ) then
        if ( codetype == do_statement ) then
            call testcov_get_count( inname, lineno, 1, count   )
            call testcov_get_count( inname, lineno, 2, skipped )
            write( 21, '(i8,a,a)' ) count, ': ', line(1)
            if ( count == 0 ) then
                write( 21, '(a8,a)'  ) '==> ', ': (do-loop never ran)'
            endif
            if ( skipped == 0 ) then
                write( 21, '(a8,a)'  ) '==> ', &
                    ': (do-loop never skipped, that is: always at least one iteration)'
            endif
            if ( nolines > 1 ) then
                write( 21, '(10x,a)'  ) ( line(i), i = 2,nolines )
            endif
        else
            write( 21, '(10x,a)'  ) ( line(i), i = 1,nolines )
        endif
    else
        if ( codetype == do_statement ) then
            do_level = do_level + 1
            do_start(do_level) = lineno
            call write_line( (/ ' ' /), lineno, 1, inname, 2 )
            call write_line( line, lineno, nolines, inname, 1 )
        else
            !
            ! TODO: should be written with lineno belonging to DO
            !
            call write_line( line, do_start(do_level), nolines, inname, 3 )
            do_level = do_level - 1

            if ( do_level < 0 ) then
                write(*,*) 'Error in handle_do: more enddo than dos!'
            endif
        endif
    endif
end subroutine handle_do

! handle_if --
!     Output for the control statements of an if-block/statement
!
! Arguments:
!     line           Line of code to instrument/report
!     lineno         Line number
!     nolines        Number of lines
!     codetype       Class to which it belongs
!
! Note:
!     No provision yet for an if-statement
!
subroutine handle_if( line, lineno, nolines, codetype )
    character(len=*), dimension(:), intent(in)    :: line
    integer, intent(in)                           :: lineno
    integer, intent(in)                           :: nolines
    integer, intent(in)                           :: codetype

    logical, dimension(100), save                 :: has_else = .false.
    integer, save                                 :: if_level = 0
    integer                                       :: count
    integer                                       :: i

    if ( report ) then
        select case (codetype)
            case ( if_statement )
                if_level = if_level + 1
                has_else(if_level) = .false.
                call testcov_get_count( inname, lineno, 1, count )
                write( 21, '(i8,a,a)' ) count, ': ', line(1)
                if ( count == 0 ) then
                    write( 21, '(a8,a)'  ) '==> ', ': (if/elseif condition was never true)'
                endif
                if ( nolines > 1 ) then
                    write( 21, '(10x,a)'  ) ( line(i), i = 2,nolines )
                endif

            case ( else_statement )
                has_else(if_level) = .true.
                call testcov_get_count( inname, lineno, 1, count )
                write( 21, '(i8,a,a)' ) count, ': ', line(1)
                if ( count == 0 ) then
                    write( 21, '(a8,a)'  ) '==> ', ': (else condition was never true)'
                endif
                if ( nolines > 1 ) then
                    write( 21, '(10x,a)'  ) ( line(i), i = 2,nolines )
                endif

            case ( endif_statement )
                write( 21, '(10x,a)'  ) ( line(i), i = 1,nolines )
                if ( .not. has_else(if_level) ) then
                    call testcov_get_count( inname, lineno, 2, count )
                    if ( count > 0 ) then
                        write( 21, '(i8,a,a)' ) count, ': (if condition false)'
                    else
                        write( 21, '(a8,a)' )   '==> ', ': (if condition always true ==> implicit else did not occur)'
                    endif
                endif
                if_level = if_level - 1

        end select
    else
        !
        ! TODO: if-statements instead of blocks
        !
        select case (codetype)
            case ( if_statement )
                if_level = if_level + 1
                has_else(if_level) = .false.
                call write_line( line, lineno, nolines, inname, 1 )

            case ( else_statement )
                has_else(if_level) = .true.
                call write_line( line, lineno, nolines, inname, 1 )

            case ( endif_statement )
                if ( .not. has_else(if_level) ) then
                    write( 21, '(a)' ) '        else'
                    call write_line( (/ ' ' /), lineno, 1, inname, 3 )
                endif
                call write_line( line, lineno, nolines, inname, 1 )
                if_level = if_level - 1

                if ( if_level < 0 ) then
                    write(*,*) 'Error in handle_if: more endifs than ifs!'
                endif

            case default
                write(*,*) 'Programming error in handle_if - unknown case'
        endselect
    endif
end subroutine handle_if

! handle_case --
!     Output for the control statements of a select-case block
!
! Arguments:
!     line           Line of code to instrument/report
!     lineno         Line number
!     nolines        Number of lines
!     codetype       Class to which it belongs
!
subroutine handle_case( line, lineno, nolines, codetype )
    character(len=*), dimension(:), intent(in)    :: line
    integer, intent(in)                           :: lineno
    integer, intent(in)                           :: nolines
    integer, intent(in)                           :: codetype

    logical, dimension(100), save                 :: has_default = .false.
    integer, save                                 :: select_level = 0
    integer                                       :: i
    integer                                       :: count

    if ( report ) then
        call testcov_get_count( inname, lineno, 1, count )
        select case (codetype)
            case ( select_statement )
                select_level = select_level + 1
                write( 21, '(10x,a)'  ) ( line(i), i = 1,nolines )

            case ( case_statement )
                write( 21, '(i8,a,a)' ) count, ': ', line(1)
                if ( count == 0 ) then
                    write( 21, '(a8,a)'  ) '==> ', ': (case never occurred)'
                endif
                if ( nolines > 1 ) then
                    write( 21, '(10x,a)'  ) ( line(i), i = 2,nolines )
                endif

            case ( casedefault_statement )
                has_default(select_level) = .true.
                write( 21, '(i8,a,a)' ) count, ': ', line(1)
                if ( count == 0 ) then
                    write( 21, '(a8,a)'  ) '==> ', ': (default case never occurred)'
                endif
                if ( nolines > 1 ) then
                    write( 21, '(10x,a)'  ) ( line(i), i = 2,nolines )
                endif

            case ( endselect_statement )
                call testcov_get_count( inname, lineno, 2, count )
                write( 21, '(i8,a,a)'  ) count, ': ', line(1)
                if ( nolines > 1 ) then
                    write( 21, '(10x,a)'  ) ( line(i), i = 2,nolines )
                endif
                if ( .not. has_default(select_level) ) then
                    if ( count == 0 ) then
                        write( 21, '(a8,a)'  ) '==> ', ': (implicit default case never occurred)'
                    endif
                endif
                select_level = select_level - 1

                if ( select_level < 0 ) then
                    write(21,*) 'Error in handle_case: more endselects than selects!'
                    write(*,*) 'Error in handle_case: more endselects than selects!'
                endif

            case default
                write(*,*) 'Programming error in handle_select - unknown case'
        endselect
    else
        !
        ! TODO: if-statements instead of blocks
        !
        select case (codetype)
            case ( select_statement )
                select_level = select_level + 1
                has_default(select_level) = .false.
                write( 21, '(a)' ) (trim(line(i)) ,i=1,nolines)

            case ( case_statement )
                call write_line( line, lineno, nolines, inname, 1 )

            case ( casedefault_statement )
                has_default(select_level) = .true.
                call write_line( line, lineno, nolines, inname, 1 )

            case ( endselect_statement )
                if ( .not. has_default(select_level) ) then
                    write( 21, '(a)' ) '        case default'
                    call write_line( line, lineno, nolines, inname, 3 )
                else
                    write( 21, '(a)' ) (trim(line(i)) ,i=1,nolines)
                endif
                select_level = select_level - 1

                if ( select_level < 0 ) then
                    write(21,*) 'Error in handle_case: more endselects than selects!'
                    write(*,*) 'Error in handle_case: more endselects than selects!'
                endif

            case default
                write(*,*) 'Programming error in handle_select - unknown case'
        endselect
    endif
end subroutine handle_case

end program


text:
In the heart of downtown Davis, Street Savvy offers stylish, quality mens clothing and accessories at a savvy price. The exclusive selection includes chic, yet price conscious younger mens lines as well as classic cuts for the established man.

The store is equipped with a friendly pubstyle bar and a 42inch plasma TV featuring ESPN. Every Friday evening at 5:00pm you can enjoy complimentary cheese, crackers and wine while the expert staff assists you with your weekend plans.
  Theres a pub inside the store? How much for a pint?  Users/ArlenAbraham
Come experience a Davis original.

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


text:
Rumor has it that new restaurants and other arriving businesses are often delayed in opening because of delays in the City performing the necessary inspections.  Can anyone confirm, deny, or elaborate?

Here is a http://www.davisenterprise.com/opinion/letters/expeditecitybusinessinspections/ Davis Enterprise letter to the editor, suggesting the following addition to Councilman Dan Wolks Action Plan for Downtown:  Speed up city reviews and inspections of new retail businesses and restaurants to get them operating in a timely fashion.

Building code side? Think of the inspection that comes from the city, code and business License. The food safety is done @ county level. It could be having permit inspection signed off on... Daubert


text:
Hi all. I never posted anything on a wiki until just tonight.  Im involved with Davis Bicycles! (and Davis bike culture in general), so Ill be back.

20070917 23:19:50 nbsp awesome work with regard to the GPS and KML data! Users/Graham.Freeman

20070918 00:27:42 nbsp those are awesome pics of Users/TedBuehler ted working on bikesUsers/JessicaRockwell


text:
 SUBROUTINE POLATES1(IPOPT,KGDSI,KGDSO,MI,MO,KM,IBI,LI,GI, &
                     NO,RLAT,RLON,IBO,LO,GO,IRET)
!$$$  SUBPROGRAM DOCUMENTATION BLOCK
!
! $Revision: 71314 $
!
! SUBPROGRAM:  POLATES1   INTERPOLATE SCALAR FIELDS (BICUBIC)
!   PRGMMR: IREDELL       ORG: W/NMC23       DATE: 96-04-10
!
! ABSTRACT: THIS SUBPROGRAM PERFORMS BICUBIC INTERPOLATION
!           FROM ANY GRID TO ANY GRID FOR SCALAR FIELDS.
!           BITMAPS ARE NOW ALLOWED EVEN WHEN INVALID POINTS ARE WITHIN
!           THE BICUBIC TEMPLATE PROVIDED THE MINIMUM WEIGHT IS REACHED.
!           OPTIONS ALLOW CHOICES BETWEEN STRAIGHT BICUBIC (IPOPT(1)=0)
!           AND CONSTRAINED BICUBIC (IPOPT(1)=1) WHERE THE VALUE IS
!           CONFINED WITHIN THE RANGE OF THE SURROUNDING 16 POINTS.
!           ANOTHER OPTION IS THE MINIMUM PERCENTAGE FOR MASK,
!           I.E. PERCENT VALID INPUT DATA REQUIRED TO MAKE OUTPUT DATA,
!           (IPOPT(2)) WHICH DEFAULTS TO 50 (IF IPOPT(2)=-1).
!           BILINEAR USED WITHIN ONE GRID LENGTH OF BOUNDARIES.
!           ONLY HORIZONTAL INTERPOLATION IS PERFORMED.
!           THE GRIDS ARE DEFINED BY THEIR GRID DESCRIPTION SECTIONS
!           (PASSED IN INTEGER FORM AS DECODED BY SUBPROGRAM W3FI63).
!           THE CURRENT CODE RECOGNIZES THE FOLLOWING PROJECTIONS:
!             (KGDS(1)=000) EQUIDISTANT CYLINDRICAL
!             (KGDS(1)=001) MERCATOR CYLINDRICAL
!             (KGDS(1)=003) LAMBERT CONFORMAL CONICAL
!             (KGDS(1)=004) GAUSSIAN CYLINDRICAL (SPECTRAL NATIVE)
!             (KGDS(1)=005) POLAR STEREOGRAPHIC AZIMUTHAL
!             (KGDS(1)=203) ROTATED EQUIDISTANT CYLINDRICAL (E-STAGGER)
!             (KGDS(1)=205) ROTATED EQUIDISTANT CYLINDRICAL (B-STAGGER)
!           WHERE KGDS COULD BE EITHER INPUT KGDSI OR OUTPUT KGDSO.
!           AS AN ADDED BONUS THE NUMBER OF OUTPUT GRID POINTS
!           AND THEIR LATITUDES AND LONGITUDES ARE ALSO RETURNED.
!           ON THE OTHER HAND, THE OUTPUT CAN BE A SET OF STATION POINTS
!           IF KGDSO(1)<0, IN WHICH CASE THE NUMBER OF POINTS
!           AND THEIR LATITUDES AND LONGITUDES MUST BE INPUT.
!           OUTPUT BITMAPS WILL ONLY BE CREATED WHEN THE OUTPUT GRID
!           EXTENDS OUTSIDE OF THE DOMAIN OF THE INPUT GRID.
!           THE OUTPUT FIELD IS SET TO 0 WHERE THE OUTPUT BITMAP IS OFF.
!
! PROGRAM HISTORY LOG:
!   96-04-10  IREDELL
! 1999-04-08  IREDELL  SPLIT IJKGDS INTO TWO PIECES
! 2001-06-18  IREDELL  INCLUDE MINIMUM MASK PERCENTAGE OPTION
! 2007-05-22  IREDELL  EXTRAPOLATE UP TO HALF A GRID CELL
! 2007-10-30  IREDELL  CORRECT NORTH POLE INDEXING PROBLEM,
!                      UNIFY MASKED AND NON-MASKED ALGORITHMS,
!                      AND SAVE WEIGHTS FOR PERFORMANCE.
! 2012-06-26  GAYNO    FIX OUT-OF-BOUNDS ERROR.  SEE NCEPLIBS
!                      TICKET #9.
! 2015-01-27  GAYNO    REPLACE CALLS TO GDSWIZ WITH NEW MERGED
!                      VERSION OF GDSWZD.
!
! USAGE:    CALL POLATES1(IPOPT,KGDSI,KGDSO,MI,MO,KM,IBI,LI,GI,
!    &                    NO,RLAT,RLON,IBO,LO,GO,IRET)
!
!   INPUT ARGUMENT LIST:
!     IPOPT    - INTEGER (20) INTERPOLATION OPTIONS
!                IPOPT(1)=0 FOR STRAIGHT BICUBIC;
!                IPOPT(1)=1 FOR CONSTRAINED BICUBIC WHERE VALUE IS
!                CONFINED WITHIN THE RANGE OF THE SURROUNDING 4 POINTS.
!                IPOPT(2) IS MINIMUM PERCENTAGE FOR MASK
!                (DEFAULTS TO 50 IF IPOPT(2)=-1)
!     KGDSI    - INTEGER (200) INPUT GDS PARAMETERS AS DECODED BY W3FI63
!     KGDSO    - INTEGER (200) OUTPUT GDS PARAMETERS
!                (KGDSO(1)<0 IMPLIES RANDOM STATION POINTS)
!     MI       - INTEGER SKIP NUMBER BETWEEN INPUT GRID FIELDS IF KM>1
!                OR DIMENSION OF INPUT GRID FIELDS IF KM=1
!     MO       - INTEGER SKIP NUMBER BETWEEN OUTPUT GRID FIELDS IF KM>1
!                OR DIMENSION OF OUTPUT GRID FIELDS IF KM=1
!     KM       - INTEGER NUMBER OF FIELDS TO INTERPOLATE
!     IBI      - INTEGER (KM) INPUT BITMAP FLAGS
!     LI       - LOGICAL*1 (MI,KM) INPUT BITMAPS (IF SOME IBI(K)=1)
!     GI       - REAL (MI,KM) INPUT FIELDS TO INTERPOLATE
!     NO       - INTEGER NUMBER OF OUTPUT POINTS (ONLY IF KGDSO(1)<0)
!     RLAT     - REAL (NO) OUTPUT LATITUDES IN DEGREES (IF KGDSO(1)<0)
!     RLON     - REAL (NO) OUTPUT LONGITUDES IN DEGREES (IF KGDSO(1)<0)
!
!   OUTPUT ARGUMENT LIST:
!     NO       - INTEGER NUMBER OF OUTPUT POINTS (ONLY IF KGDSO(1)>=0)
!     RLAT     - REAL (MO) OUTPUT LATITUDES IN DEGREES (IF KGDSO(1)>=0)
!     RLON     - REAL (MO) OUTPUT LONGITUDES IN DEGREES (IF KGDSO(1)>=0)
!     IBO      - INTEGER (KM) OUTPUT BITMAP FLAGS
!     LO       - LOGICAL*1 (MO,KM) OUTPUT BITMAPS (ALWAYS OUTPUT)
!     GO       - REAL (MO,KM) OUTPUT FIELDS INTERPOLATED
!     IRET     - INTEGER RETURN CODE
!                0    SUCCESSFUL INTERPOLATION
!                2    UNRECOGNIZED INPUT GRID OR NO GRID OVERLAP
!                3    UNRECOGNIZED OUTPUT GRID
!
! SUBPROGRAMS CALLED:
!   GDSWZD       GRID DESCRIPTION SECTION WIZARD
!   IJKGDS0      SET UP PARAMETERS FOR IJKGDS1
!   (IJKGDS1)    RETURN FIELD POSITION FOR A GIVEN GRID POINT
!   POLFIXS      MAKE MULTIPLE POLE SCALAR VALUES CONSISTENT
!
! ATTRIBUTES:
!   LANGUAGE: FORTRAN 90
!
!$$$
!
 USE GDSWZD_MOD
!
 IMPLICIT NONE
!
 INTEGER,                INTENT(IN   ):: IPOPT(20),KGDSI(200)
 INTEGER,                INTENT(IN   ):: KGDSO(200),MI,MO,KM
 INTEGER,                INTENT(IN   ):: IBI(KM)
 INTEGER,                INTENT(INOUT):: NO
 INTEGER,                INTENT(  OUT):: IRET, IBO(KM)
!
 LOGICAL*1,              INTENT(IN   ):: LI(MI,KM)
 LOGICAL*1,              INTENT(  OUT):: LO(MO,KM)
!
 REAL,                   INTENT(IN   ):: GI(MI,KM)
 REAL,                   INTENT(INOUT):: RLAT(MO),RLON(MO)
 REAL,                   INTENT(  OUT):: GO(MO,KM)
!
 REAL,                   PARAMETER    :: FILL=-9999.
!
 INTEGER                              :: IJKGDSA(20)
 INTEGER                              :: IJX(4),IJY(4)
 INTEGER                              :: MCON,MP,N,I,J,K
 INTEGER                              :: NK,NV,IJKGDS1
 INTEGER,            SAVE             :: KGDSIX(200)=-1,KGDSOX(200)=-1
 INTEGER,            SAVE             :: NOX=-1,IRETX=-1
 INTEGER,            ALLOCATABLE,SAVE :: NXY(:,:,:),NC(:)
!
 REAL                                 :: PMP,XIJ,YIJ,XF,YF
 REAL                                 :: G,W,GMIN,GMAX
 REAL                                 :: WX(4),WY(4)
 REAL                                 :: XPTS(MO),YPTS(MO)
 REAL,               ALLOCATABLE,SAVE :: RLATX(:),RLONX(:),WXY(:,:,:)
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  SET PARAMETERS
 IRET=0
 MCON=IPOPT(1)
 MP=IPOPT(2)
 IF(MP.EQ.-1.OR.MP.EQ.0) MP=50
 IF(MP.LT.0.OR.MP.GT.100) IRET=32
 PMP=MP*0.01
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  SAVE OR SKIP WEIGHT COMPUTATION
 IF(IRET.EQ.0.AND.(KGDSO(1).LT.0.OR.ANY(KGDSI.NE.KGDSIX).OR.ANY(KGDSO.NE.KGDSOX))) THEN
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  COMPUTE NUMBER OF OUTPUT POINTS AND THEIR LATITUDES AND LONGITUDES.
   IF(KGDSO(1).GE.0) THEN
     CALL GDSWZD(KGDSO, 0,MO,FILL,XPTS,YPTS,RLON,RLAT,NO)
     IF(NO.EQ.0) IRET=3
   ENDIF
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  LOCATE INPUT POINTS
   CALL GDSWZD(KGDSI,-1,NO,FILL,XPTS,YPTS,RLON,RLAT,NV)
   IF(IRET.EQ.0.AND.NV.EQ.0) IRET=2
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  ALLOCATE AND SAVE GRID DATA
   KGDSIX=KGDSI
   KGDSOX=KGDSO
   IF(NOX.NE.NO) THEN
     IF(NOX.GE.0) DEALLOCATE(RLATX,RLONX,NC,NXY,WXY)
     ALLOCATE(RLATX(NO),RLONX(NO),NC(NO),NXY(4,4,NO),WXY(4,4,NO))
     NOX=NO
   ENDIF
   IRETX=IRET
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  COMPUTE WEIGHTS
   IF(IRET.EQ.0) THEN
     CALL IJKGDS0(KGDSI,IJKGDSA)
!$OMP PARALLEL DO PRIVATE(N,XIJ,YIJ,IJX,IJY,XF,YF,J,I,WX,WY)
     DO N=1,NO
       RLONX(N)=RLON(N)
       RLATX(N)=RLAT(N)
       XIJ=XPTS(N)
       YIJ=YPTS(N)
       IF(XIJ.NE.FILL.AND.YIJ.NE.FILL) THEN
         IJX(1:4)=FLOOR(XIJ-1)+(/0,1,2,3/)
         IJY(1:4)=FLOOR(YIJ-1)+(/0,1,2,3/)
         XF=XIJ-IJX(2)
         YF=YIJ-IJY(2)
         DO J=1,4
           DO I=1,4
             NXY(I,J,N)=IJKGDS1(IJX(I),IJY(J),IJKGDSA)
           ENDDO
         ENDDO
         IF(MINVAL(NXY(1:4,1:4,N)).GT.0) THEN
!  BICUBIC WHERE 16-POINT STENCIL IS AVAILABLE
           NC(N)=1
           WX(1)=XF*(1-XF)*(2-XF)/(-6.)
           WX(2)=(XF+1)*(1-XF)*(2-XF)/2.
           WX(3)=(XF+1)*XF*(2-XF)/2.
           WX(4)=(XF+1)*XF*(1-XF)/(-6.)
           WY(1)=YF*(1-YF)*(2-YF)/(-6.)
           WY(2)=(YF+1)*(1-YF)*(2-YF)/2.
           WY(3)=(YF+1)*YF*(2-YF)/2.
           WY(4)=(YF+1)*YF*(1-YF)/(-6.)
         ELSE
!  BILINEAR ELSEWHERE NEAR THE EDGE OF THE GRID
           NC(N)=2
           WX(1)=0
           WX(2)=(1-XF)
           WX(3)=XF
           WX(4)=0
           WY(1)=0
           WY(2)=(1-YF)
           WY(3)=YF
           WY(4)=0
         ENDIF
         DO J=1,4
           DO I=1,4
             WXY(I,J,N)=WX(I)*WY(J)
           ENDDO
         ENDDO
       ELSE
         NC(N)=0
       ENDIF
     ENDDO
   ENDIF
 ENDIF
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  INTERPOLATE OVER ALL FIELDS
 IF(IRET.EQ.0.AND.IRETX.EQ.0) THEN
   IF(KGDSO(1).GE.0) THEN
     NO=NOX
     DO N=1,NO
       RLON(N)=RLONX(N)
       RLAT(N)=RLATX(N)
     ENDDO
   ENDIF
!$OMP PARALLEL DO PRIVATE(NK,K,N,G,W,GMIN,GMAX,J,I)
   DO NK=1,NO*KM
     K=(NK-1)/NO+1
     N=NK-NO*(K-1)
     IF(NC(N).GT.0) THEN
       G=0
       W=0
       IF(MCON.GT.0) GMIN=HUGE(GMIN)
       IF(MCON.GT.0) GMAX=-HUGE(GMAX)
       DO J=NC(N),5-NC(N)
         DO I=NC(N),5-NC(N)
           IF(NXY(I,J,N).GT.0)THEN
             IF(IBI(K).EQ.0.OR.LI(NXY(I,J,N),K))THEN
               G=G+WXY(I,J,N)*GI(NXY(I,J,N),K)
               W=W+WXY(I,J,N)
               IF(MCON.GT.0) GMIN=MIN(GMIN,GI(NXY(I,J,N),K))
               IF(MCON.GT.0) GMAX=MAX(GMAX,GI(NXY(I,J,N),K))
             ENDIF
           ENDIF
         ENDDO
       ENDDO
       LO(N,K)=W.GE.PMP
       IF(LO(N,K)) THEN
         GO(N,K)=G/W
         IF(MCON.GT.0) GO(N,K)=MIN(MAX(GO(N,K),GMIN),GMAX)
       ELSE
         GO(N,K)=0.
       ENDIF
     ELSE
       LO(N,K)=.FALSE.
       GO(N,K)=0.
     ENDIF
   ENDDO
   DO K=1,KM
     IBO(K)=IBI(K)
     IF(.NOT.ALL(LO(1:NO,K))) IBO(K)=1
   ENDDO
   IF(KGDSO(1).EQ.0) CALL POLFIXS(NO,MO,KM,RLAT,RLON,IBO,LO,GO)
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 ELSE
   IF(IRET.EQ.0) IRET=IRETX
   IF(KGDSO(1).GE.0) NO=0
 ENDIF
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 END SUBROUTINE POLATES1

max_stars_repo_name: AerodyneLabs/Stratocast

max_issues_repo_name: AerodyneLabs/Stratocast

max_forks_repo_name: AerodyneLabs/Stratocast


!
!     ---------------------------------------------------------------
!     M E S
!     ---------------------------------------------------------------
!
!     Written by G. Gaigalas,                                      *
!     Vilnius,  Lithuania                          December 1993   *
!
      SUBROUTINE MES(I)
!...Translated by Pacific-Sierra Research 77to90  4.3E  10:52:55  11/16/01
!...Switches:
      IMPLICIT NONE
!-----------------------------------------------
!   D u m m y   A r g u m e n t s
!-----------------------------------------------
      INTEGER , INTENT(IN) :: I
!-----------------------------------------------
!   L o c a l   V a r i a b l e s
!-----------------------------------------------
      INTEGER :: J
      CHARACTER, DIMENSION(6) :: STRING5*10
!-----------------------------------------------
      DATA STRING5/ '   I T L S', ' I T L S 2', ' I T L S 3', 'AW P 1 L S', &
         'WA P 1 L S', '       W 1'/
      IF (I > 50) THEN
         J = I - 50
         WRITE (6, '(A)') ' error in func./sub. '
         WRITE (6, '(20X,A10)') STRING5(J)
         WRITE (6, '(A)') ' susimaise f sluoksnio termu kodavimas  '
      ELSE
         WRITE (6, '(A)') ' yra daugiau nei 2 ele. sluoks. f,g,h,i,k,l,m'
         WRITE (6, '(3X,I5)') I
         SELECT CASE (I)
         CASE (1)
            WRITE (6, '(A)') ' error in Subroutine   W 1 G '
         CASE (2)
            WRITE (6, '(A)') ' error in Subroutine   W 1 '
         CASE (11)
            WRITE (6, '(A)') ' error in Function     I T L S  '
         CASE (12)
            WRITE (6, '(A)') ' error in Function     I T L S 2  '
         CASE (13)
            WRITE (6, '(A)') ' error in Function     I T L S 3  '
         CASE (30)
            WRITE (6, '(A)') ' error in Subroutine   A 1 A 2 A 3 A 4 L S '
         CASE (31)
            WRITE (6, '(A)') ' error in Subroutine   A 1 A 2 L S '
         CASE (32)
            WRITE (6, '(A)') ' error in Subroutine   A 1 A 2 W 3 L S '
         CASE (33)
            WRITE (6, '(A)') ' error in Subroutine   A 1 A W 2 L S '
         CASE (34)
            WRITE (6, '(A)') ' error in Subroutine   W A 1 A 2 L S '
         CASE (35)
            WRITE (6, '(A)') ' error in Subroutine   W 1 W 2 L S '
         CASE DEFAULT
            WRITE (6, '(A)') ' error in unknown Subroutine  '
         END SELECT
      ENDIF
      WRITE (6, '(A)') ' Contact to   G. Gediminas please ! ! ! '
      STOP
      END SUBROUTINE MES


text: #redirect EP Copier Cooperative de Alley


SUBROUTINE empirical_cor (Xmatrix, Bias_0, CPR_0, Bias_apos, CPR_apos)

! ----------------------------------------------------------------------
! SUBROUTINE: empirical_cor.f03
! ----------------------------------------------------------------------
! Purpose:
!  Orbit empirical parameters initialization
! ----------------------------------------------------------------------
! Input arguments:
! - Bias_0: 	Bias parameters acceleration apriori values
! - CPR_0: 		Cycle per revolution parameters (C and S coefficients) apriori values
!
! Output arguments:
! - Bias_apos: 	Bias parameters acceleration aposteriori values
! - CPR_apos: 	Cycle per revolution parameters (C and S coefficients) aposteriori values
! ----------------------------------------------------------------------
! Author :	Dr. Thomas Papanikolaou
!			Geoscience Australia, CRC-SI
! Created:	30 November 2018
! ----------------------------------------------------------------------

      USE mdl_precision
      USE mdl_num
      USE mdl_param
      use pod_yaml
      IMPLICIT NONE

! ----------------------------------------------------------------------
! Dummy arguments declaration
! ----------------------------------------------------------------------
! IN
      REAL (KIND = prec_d), DIMENSION(:,:), ALLOCATABLE, INTENT(IN) :: Xmatrix
	  REAL (KIND = prec_q), INTENT(IN) :: Bias_0(3)
	  REAL (KIND = prec_q), INTENT(IN) :: CPR_0(3,2)
! ----------------------------------------------------------------------
! OUT
	  REAL (KIND = prec_q), INTENT(OUT) :: Bias_apos(3)
	  REAL (KIND = prec_q), INTENT(OUT) :: CPR_apos(3,2)
! ----------------------------------------------------------------------

! ----------------------------------------------------------------------
! Local variables declaration
! ----------------------------------------------------------------------
      REAL (KIND = prec_d) :: Bias_corr(3), CPR_corr(3,2)
! ----------------------------------------------------------------------

! ----------------------------------------------------------------------
! FIXME: variable initialisation
Bias_corr = 0.d0
CPR_corr = 0.d0

! ----------------------------------------------------------------------
If (BTEST(yml_srp_parameters, EMP_R_Bias - one) .and. BTEST(yml_srp_parameters, EMP_R_CPR - one)) Then
! Bias & CPR
Bias_corr = Xmatrix(7:9,1)
CPR_corr(1,1) = Xmatrix(10,1)
CPR_corr(1,2) = Xmatrix(11,1)
CPR_corr(2,1) = Xmatrix(12,1)
CPR_corr(2,2) = Xmatrix(13,1)
CPR_corr(3,1) = Xmatrix(14,1)
CPR_corr(3,2) = Xmatrix(15,1)

else If (.not. BTEST(yml_srp_parameters, EMP_R_Bias - one) .and. BTEST(yml_srp_parameters, EMP_R_CPR - one)) Then
! CPR
CPR_corr(1,1) = Xmatrix(7,1)
CPR_corr(1,2) = Xmatrix(8,1)
CPR_corr(2,1) = Xmatrix(9,1)
CPR_corr(2,2) = Xmatrix(10,1)
CPR_corr(3,1) = Xmatrix(11,1)
CPR_corr(3,2) = Xmatrix(12,1)

else If (BTEST(yml_srp_parameters, EMP_R_Bias - one) .and. .not. BTEST(yml_srp_parameters, EMP_R_CPR - one)) Then
! Bias
Bias_corr = Xmatrix(7:9,1)

End If
! ----------------------------------------------------------------------

Bias_apos = Bias_0 + Bias_corr
!print *, "Bias_accel_aposteriori", Bias_accel_aposteriori
!print *, "Bias_accel_glb", Bias_accel_glb
!print *, "Bias_corr", Bias_corr

CPR_apos = CPR_0 + CPR_corr
!print *, "CPR_CS_aposteriori", CPR_CS_aposteriori
!print *, "CPR_CS_glb", CPR_CS_glb
!print *, "CPR_corr", CPR_corr

END SUBROUTINE


text:
The Alamo West Centre is a series of two story office buildings in between D Street and C Street while facing  2nd Street.

 First Floor Offices
 Davisville Travel
 Fidelity National Title
 First Northern Bank
 Mortgage Resource Network Inc.
 There might be more

 Second Floor Offices
 Cedant Inc.
 Ceres Corporation
 Davis Shambhala Meditation Center
 Institute of Behavior Assessment
 Law Office of Stephen M. Reynolds
 There might be more


!
!     -------------------------------------------------------------
!      S S 1 2 1 2
!     -------------------------------------------------------------
!                                                                  *
!                                                                  *
!     THIS PACKAGE DETERMINES THE VALUES OF MATRIX ELEMENTS        *
!     OF ONE PARTICLE OPERATOR IN CASE :           N'1 = N1        *
!                                                                  *
!                                                                  *
!                                                                  *
!     Written by G. Gaigalas,                                      *
!     Vanderbilt University,  Nashville             October 1996   *
!
!
!GG elektrostatine
!      CALL COULOMBLS(LJ(IA),LJ(IB),LJ(IA),LJ(IB),KL1,A1)
!      A1=A1*TWO*DSQRT(DBLE(2*KL1+1))
!GG elektrostatine
      SUBROUTINE SS1212(IG, KL, IA, IB)
!-----------------------------------------------
!   M o d u l e s
!-----------------------------------------------
      USE vast_kind_param, ONLY:  DOUBLE
      USE CONSTS_C
      USE MEDEFN_C
      USE DIAGNL_C
!...Translated by Pacific-Sierra Research 77to90  4.3E  12:25:04  11/16/01
!...Switches:
!-----------------------------------------------
!   I n t e r f a c e   B l o c k s
!-----------------------------------------------
      USE ss1_I
      USE two12_I
      USE savenon_I
      IMPLICIT NONE
!-----------------------------------------------
!   D u m m y   A r g u m e n t s
!-----------------------------------------------
      INTEGER  :: IG
      INTEGER  :: KL
      INTEGER  :: IA
      INTEGER  :: IB
!-----------------------------------------------
!   L o c a l   V a r i a b l e s
!-----------------------------------------------
      INTEGER :: KLP, LA, LB
      REAL(DOUBLE) :: A1, C
!-----------------------------------------------
      KLP = KL + 2
!      IF(IG-1.LT.KLP) RETURN
      LA = IJFUL(IA)
      LB = IJFUL(IB)
      CALL SS1 (LJ(IA), LJ(IB), KLP, KL, A1)
      IF (DABS(A1) > EPS) THEN
         CALL TWO12 (KLP, 1, KL, 1, 2, IA, IB, C)
         C = C*A1
         IF (DABS(C) > EPS) CALL SAVENON (8, C, KL, LA, LB, LA, LB, JA, JB, 0)
      ENDIF
      CALL SS1 (LJ(IB), LJ(IA), KLP, KL, A1)
      IF (DABS(A1) <= EPS) RETURN
      CALL TWO12 (KL, 1, KLP, 1, 2, IA, IB, C)
      C = C*A1
      IF (DABS(C) > EPS) CALL SAVENON (8, C, KL, LB, LA, LB, LA, JA, JB, 0)
      RETURN
      END SUBROUTINE SS1212


text:
C    @(#)kmpqpn.f	20.3 2/13/96
        function kmpqpn(m,n)
C
      include 'ipfinc/parametr.inc'

      include 'ipfinc/komps.inc'
      include 'ipfinc/optim1.inc'
C
        kmpqpn=kqpent(1,m)-kqpent(1,n)
        if(kmpqpn.eq.0 .and. m.ne.n) ndup=1
        return
        end


text:
*----------------------------------------------------------------------*
      integer function msgmdid(iocc,msd,gmd,nsym)
*----------------------------------------------------------------------*
*     calculate integer-valued ID of (ms,Gamma) distribution
*----------------------------------------------------------------------*
      implicit none

      include 'opdim.h'
      include 'hpvxseq.h'

      integer, intent(in) ::
     &     iocc(ngastp,2), msd(ngastp,2), gmd(ngastp,2), nsym

      integer ::
     &     ipatt, ibase, ica, ihpvdx, ihpv

      ipatt = 0
      ibase = 1
      do ica = 1,2
        do ihpvdx = 1, ngastp
          ihpv = hpvxseq(ihpvdx)
          if (iocc(ihpv,ica).eq.0) cycle
          ipatt = ipatt +
     &        (((iocc(ihpv,ica)-msd(ihpv,ica))/2)*nsym+
     &                                  gmd(ihpv,ica)-1)*ibase
          ibase = ibase*(iocc(ihpv,ica)+1)*nsym
        end do
      end do

      msgmdid = ipatt

      return
      end

max_stars_repo_name: ak-ustutt/GeCCo-public

max_issues_repo_name: ak-ustutt/GeCCo-public

max_forks_repo_name: ak-ustutt/GeCCo-public


text:
C Copyright(C) 1999-2020 National Technology & Engineering Solutions
C of Sandia, LLC (NTESS).  Under the terms of Contract DE-NA0003525 with
C NTESS, the U.S. Government retains certain rights in this software.
C
C See packages/seacas/LICENSE for details

C=======================================================================
      SUBROUTINE MPRESE()
      COMMON /MAP/MODEL(4,4),VIEW(4,4),PROJ(4,4),CPNEAR,CPFAR,VWPORT(4),
     *       MVP(4,4),VP(4,4),CPLINE(2,2,10),CPPLAN(2,3,10),PEYE(3),
     *       PLOOK(3),ETWIST,NCPLIN,NCPLAN,TMAT1(4,4),TMAT2(4,4),
     *       TMAT3(4,4),TVEC1(4),TVEC2(4),TVEC3(4),TVEC4(4),TARR1(32),
     *       TARR2(32),TARR3(32),TARR4(32),TARR5(32),TARR6(32),
     *       TARR7(32),TARR8(32)
      REAL MODEL,MVP

      CALL MXIDEN(4,MODEL)
      CALL MXMULT(4,MODEL,VP,MVP)
      RETURN

      END

max_stars_repo_name: jschueller/seacas

max_issues_repo_name: jschueller/seacas

max_forks_repo_name: jschueller/seacas


text:
Dylan Beaudette is a Ph.D Candidate in the Soils and Biogeochemistry Graduate Group, in the College of Ag. here at UCD. Current research topics include, but are not limited to:
  soillandscape relationship modeling
  quantitative soil property modeling
  http://casoilresource.lawr.ucdavis.edu/drupal/node/27 an online soil survey
  Soil Taxonomy related studies
  http://casoilresource.lawr.ucdavis.edu/drupal/node/95 utilization of open source software in the sciences
  mapping soils at Pinnacles National Monument
  general soil survey related matters

Some other things that tend to fill Dylans time (when he sould be doing other things) include, but are not limited to:
  programming in SQL/PHP/Python/Bash/AWK
  PHP/GD experiments
  http://trac.osgeo.org/grass/ GRASS source code hacking
  reading pages on the Davis WIki and other such online time sinks
  electronics / microcontroller programming
  cartography and map making
  pointless projects Davis isnt flat 1 WiFi Hotspots 2
  outdoor adventuring in the central Sierra
  Hobie Cat sailing

Some things that Dylan would like to be doing:
  APRS and other HAM radio things (any other hams in Davis want to start a club?)
  working on his research
  soil science education and outreach

Your sign map is an awesome idea!  You need to get on it this weekend and map out further!

Thanks! Dylan

Dylan, a few months ago you posted on the lugod mailing list asking if people would be interested in a GRASS users group. Did anything ever come of this? I would be interested. Also, I have some questions about getting California maps that maybe you can answer. Email me at nschmalenberger at fastmail.fm . Thanks. Users/NickSchmalenberger

20051106 18:48:32 nbsp The sign map is cool! Did you do it using GPS, or did you keep a text list, or annotate a map you carried, or ? (Now, could GRASS help with it?) Users/SteveDavison

20051123 09:20:03 nbsp Thanks. actually I collected the data with a little GPS, while riding my bike. I kept track of yes/no in a note pad. Then I dumped the GPS locations to a text file and added the yes/no data. GRASS was used to convert this text file into a vector dataset that I could then overlay with other things like the aerial photo and bike trail. Users/DylanBeaudette  NICE!  Users/rocksanddirt

20071112 17:11:26 nbsp That elevation map was cool.  I always figured Davis had a slight slope to it.  Users/EdWins

20080305 10:16:36 nbsp Hi Dylan, Im just responding to your comment on the La Buena Vida page. Whether HOAs serve a purpose or not in individually owned neighborhoods is a separate issue, but how would you suggest common property be managed in a development such as this? Users/MattCzarnowski


text:
Nate is Canadian.

Nate is loud and boisterous.

Nate owns the Lukes Lounge in Dixon.


text:
program test_splines
use types, only: dp
use constants, only: pi
use utils, only: assert
use mesh, only: meshexp, linspace
use splines, only: spline3, iix, iixmin, spline3ders
implicit none

real(dp), parameter :: eps = 1e-12
real(dp), allocatable :: x(:), y(:), x2(:), y2(:)
integer :: N, N2, i, ip
real(dp) :: t1, t2

! Basics:
allocate(x(5), y(5))
x = [0._dp, 1._dp, 2._dp, 3._dp, 4._dp]
y = [0._dp, 1._dp, 2._dp, 3._dp, 4._dp]
call assert(all(abs(spline3(x, y, [0._dp]) - [0._dp]) < eps))
call assert(all(abs(spline3(x, y, [0._dp, 0.5_dp, 1._dp, 10._dp]) - &
        [0._dp, 0.5_dp, 1._dp, 10._dp]) < eps))
call assert(all(abs(spline3(x, y, [10._dp, 1._dp, 0._dp, 0.5_dp, 1._dp, &
        10._dp]) - [10._dp, 1._dp, 0._dp, 0.5_dp, 1._dp, 10._dp]) < eps))

x = linspace(0.0_dp, pi, 5)
y = sin(x)
allocate(x2(100))
x2 = linspace(0.0_dp, pi, 100)
call assert(all(abs(spline3(x, y, x2) - sin(x2)) < 1.2e-2_dp))
deallocate(x, y, x2)

! Test iix correctness:
N = 100000
N2 = 1000
allocate(x(N+1), x2(N2+1))
x = meshexp(1e-7_dp, 50._dp, 1e5_dp, N)
x2 = meshexp(0._dp, 51._dp, 200._dp, N2)
call cpu_time(t1)
do i = 1, N2+1
    ip = iix(x2(i), x)
    if (x2(i) < x(1)) then
        call assert(ip == 1)
    elseif (x2(i) > x(size(x))) then
        call assert(ip == size(x)-1)
    else
        call assert(x(ip) <= x2(i) .and. x2(i) <= x(ip+1))
    end if
end do
call cpu_time(t2)
print *, "Pure iix:", t2-t1

! Test faster indexing:
call cpu_time(t1)
ip = 0
do i = 1, N2+1
    ip = iixmin(x2(i), x, ip)
    if (x2(i) < x(1)) then
        call assert(ip == 1)
    elseif (x2(i) > x(size(x))) then
        call assert(ip == size(x)-1)
    else
        call assert(x(ip) <= x2(i) .and. x2(i) <= x(ip+1))
    end if
end do
call cpu_time(t2)
print *, "Reducing size iix:", t2-t1
deallocate(x, x2)

! Correctness. This tests that BC at each end is correct:
N = 10
N2 = 1000
allocate(x(N+1), y(N+1), x2(N2+1), y2(N2+1))
x = meshexp(0._dp, pi/2, 1._dp, N)
y = sin(x)
x2 = meshexp(-pi/2, 3*pi/4, 1._dp, N2)
y2 = spline3(x, y, x2)
call assert(abs(y2(1) - sin(x2(1))) < 5e-2_dp)
call assert(abs(y2(size(y2)) - sin(x2(size(x2)))) < 5e-2_dp)
! 1st derivative
call spline3ders(x, y, x2, dynew=y2)
call assert(abs(y2(1) - cos(x2(1))) < 0.2_dp)
call assert(abs(y2(size(y2)) - cos(x2(size(x2)))) < 0.2_dp)
! 2nd derivative
call spline3ders(x, y, x2, d2ynew=y2)
call assert(abs(y2(1) + sin(x2(1))) < 0.52_dp)
call assert(abs(y2(size(y2)) + sin(x2(size(x2)))) < 0.52_dp)

x = meshexp(0._dp, 1._dp, 1._dp, N)
y = x**2
x2 = meshexp(-10._dp, 10._dp, 1._dp, N2)
y2 = spline3(x, y, x2)
call assert(abs(y2(1) - x2(1)**2) < 1e-10_dp)
call assert(abs(y2(size(y2)) - x2(size(x2))**2) < 1e-10_dp)
! 1st derivative
call spline3ders(x, y, x2, dynew=y2)
call assert(abs(y2(1) - 2*x2(1)) < 1e-10_dp)
call assert(abs(y2(size(y2)) - 2*x2(size(x2))) < 1e-10_dp)
! 2nd derivative
call spline3ders(x, y, x2, d2ynew=y2)
call assert(maxval(abs(y2 - 2)) < 1e-10_dp)
deallocate(x, y, x2, y2)

! Test speed
N  = 1000
N2 = 10000
allocate(x(N+1), y(N+1), x2(N2+1), y2(N2+1))
x = meshexp(1e-7_dp, 50._dp, 1e5_dp, N)
x2 = meshexp(0._dp, 51._dp, 200._dp, N2)
y = sin(x)
call cpu_time(t1)
y2 = spline3(x, y, x2)
call cpu_time(t2)
print *, "Spline time:", t2-t1

end program

max_stars_repo_name: wave-scattering/fortran-utils


text:
c-----------------------------------------------------------------------
c
c-----------------------------------------------------------------------
c
c                           LIMEX4_2_Dense
c
c-----------------------------------------------------------------------
c
c     This is a collection of subroutines for LIMEX, the extrapolation
c     integrator for  the solution  of linearly-implicit differential-
c     algebraic systems of the form
c
c          B (t,y) * y' (t) = f (t,y)
c
c     with B a (n,n)-matrix of rank less or equal n.
c
c     This  collection contains  routines for  the computation and the
c     evaluation of right derivatives and Hermite polynomials.
c
c-----------------------------------------------------------------------
c
c     Copyright (C) 2000, Konrad-Zuse-Zentrum fuer Informationstechnik
c     Berlin (ZIB)
c     ALL RIGHTS RESERVED
c
c     Written by:
c
c     R. Ehrig, U. Nowak,
c     Konrad-Zuse-Zentrum fuer Informationstechnik Berlin (ZIB)
c     Takustrasse 7
c     D-14195 Berlin-Dahlem
c
c     phone : +49-30-84185-0
c     fax   : +49-30-84185-125
c     e-mail: ehrig@zib.de, nowak@zib.de
c     URL   : http://www.zib.de
c
c-----------------------------------------------------------------------
c
c            **************************************************
c            **                                              **
c            **    This is version 4.2 of March, 17, 2000    **
c            **                                              **
c            **************************************************
c
c-----------------------------------------------------------------------
c
c     Overview of current versions:
c
c     4.2A1  Non-sparse dense  or  banded  Jacobians. Direct  solvers:
c            LAPACK routines. Based on BLAS and LAPACK routines.
c
c     4.2A2  Non-sparse dense or banded Jacobians. Direct solvers: NAG
c            routines. Based on NAG routines (Mark 16).
c
c     4.2B1  Sparse Jacobians. Direct  solvers: MA28 routines from the
c            Harwell subroutine library, iterative  solvers: GMRES and
c            BICGSTAB  with a  variable ILU  preconditioner. Based  on
c            BLAS and LAPACK routines.
c
c     4.2B2  Sparse Jacobians. Direct solvers: NAG routines, iterative
c            solvers: GMRES  and  BICGSTAB  with a  variable ILU  pre-
c            conditioner. Based on NAG routines (Mark 16).
c
c     Versions with other solver routines are available on request.
c
c-----------------------------------------------------------------------
c
c     NOTICE: "The LIMEX  program may be used SOLELY  for educational,
c     research, and benchmarking  purposes by no-profit organizations.
c     Commercial and other organizations  may make use of LIMEX SOLELY
c     for benchmarking  purposes only. LIMEX may  be modified by or on
c     behalf of the user for  such use but  at no time  shall LIMEX or
c     any such  modified version  of LIMEX become the  property of the
c     user. LIMEX  is provided  without warranty  of any  kind, either
c     expressed  or implied.  Neither the  authors nor their employers
c     shall be liable  for any direct or consequential  loss or damage
c     whatsoever  arising out  of the  use or  misuse of  LIMEX by the
c     user. LIMEX  must not  be sold. You  may make  copies  of LIMEX,
c     but  this NOTICE  and the  Copyright notice  must appear  in all
c     copies. Any  other  use  of LIMEX  requires written  permission.
c     Your use of LIMEX is an implicit agreement to these conditions."
c
c     LIMEX is even  available for commercial use. Contact the authors
c     who will provide details of price and condition of use.
c
c-----------------------------------------------------------------------
c
c     Contents:
c
c     Comp_Herm : computes right derivatives and the coefficient of an
c                 Hermite interpolation polynomial.
c
c     Eval_Herm : evaluates an Hermite polynomial.
c
c-----------------------------------------------------------------------
c
      subroutine Comp_Herm ( n, Max_Size, Dense, k, ipt, nj, Work )
c
      implicit none
c
c-----------------------------------------------------------------------
c
c     This  subroutine computes  approximations  for right derivatives
c     the coefficients of a Hermite interpolation polynomial.
c
c-----------------------------------------------------------------------
c
c     Copyright (C) 2000, Konrad-Zuse-Zentrum fuer Informationstechnik
c     Berlin (ZIB)
c     ALL RIGHTS RESERVED
c
c     Written by:
c
c     R. Ehrig, U. Nowak,
c     Konrad-Zuse-Zentrum fuer Informationstechnik Berlin (ZIB)
c     Takustrasse 7
c     D-14195 Berlin-Dahlem
c
c     phone : +49-30-84185-0
c     fax   : +49-30-84185-125
c     e-mail: ehrig@zib.de, nowak@zib.de
c     URL   : http://www.zib.de
c
c-----------------------------------------------------------------------
c
c     Arguments:
c
c-----------------------------------------------------------------------
c
c     n          Integer variable, must be set by caller on input, not
c                modified. Size of the system.
c
c     Max_Size   Integer variable, must be set by caller on input, not
c                modified. Leading  dimension of  the  array  Dense as
c                declared in LIMEX.
c
c     Dense      Real array of size at least Max_Size * ipt(k+1), must
c                be set by caller on input. The values of intermediate
c                solutions resp. their differences  as computed in the
c                current integration step. On exit the coefficients of
c                the Hermite interpolation poynomial.
c
c     k          Integer variable, must be set by caller on input. The
c                currrent order of the extrapolation scheme.
c
c     ipt        Integer  array of size at least k + 1, must be set by
c                caller on  input. Pointer array  to the  intermediate
c                solutions.
c
c     nj         Integer array  of size at least k + 1, must be set by
c                caller on input. The stepsize dividing sequence.
c
c     Work       Real  array of  size at  least n * (k+1), need not be
c                ste by caller on input. Work array.
c
c-----------------------------------------------------------------------
c
c     Subroutines and functions called
c
c-----------------------------------------------------------------------
c
c     From the BLAS library:
c     ----------------------
c
c     daxpy      Adds a scalar multiple of a vector to another vector
c
c     dcopy      Copies real vectors
c
c-----------------------------------------------------------------------
c
c     Declaration of the arguments
c
c-----------------------------------------------------------------------
c
      integer            n, Max_Size, k, ipt ( * ), nj ( * )
c
      double precision   Dense ( Max_Size, * ), Work ( n, * )
c
c-----------------------------------------------------------------------
c
c     Local scalars
c
c-----------------------------------------------------------------------
c
      integer            i, ifac, j, j1, j2
c
      double precision   fac, tmp
c
c-----------------------------------------------------------------------
c
c     Internal parameters
c
c-----------------------------------------------------------------------
c
      double precision   one
c
      parameter ( one = 1.0d0 )
c
c-----------------------------------------------------------------------
c
c     Main loop until the current order k.
c
c-----------------------------------------------------------------------
c
      ifac = 1
c
      do j = 1, k
c
         ifac = ifac * j
         fac  = one / dble ( ifac )
c
c-----------------------------------------------------------------------
c
c     Store the  normalized right  differences in the  first column of
c     the extrapolation tableau.
c
c-----------------------------------------------------------------------
c
         do j1 = j, k + 1
c
            tmp = fac * dble ( nj(j1) ** j )
            j2  = ipt(j1)
c
            do i = 1, n
               Work(i,j1) = tmp * Dense(i,j2)
            end do
c
         end do
c
c-----------------------------------------------------------------------
c
c     Extrapolation of the right differences.
c
c-----------------------------------------------------------------------
c
         do j1 = j, k
c
            do j2 = j1, j, - 1
c
               tmp = one / ( dble ( nj(j1+1) ) / dble ( nj(j2) ) - one )
c
               do i = 1, n
                  Work(i,j2) =   Work(i,j2+1)
     2                         + tmp * ( Work(i,j2+1) - Work(i,j2) )
               end do
c
            end do
c
         end do
c
c-----------------------------------------------------------------------
c
c     Compute higher order differences.
c
c-----------------------------------------------------------------------
c
         if ( j .lt. k ) then
c
            do j1 = j, k
               do i = ipt(j1+1), ipt(j1+1) - j1 + j, - 1
                  call daxpy ( n, -one, Dense(1,i-1), 1, Dense(1,i),
     2                         1 )
               end do
            end do
c
         end if
c
c-----------------------------------------------------------------------
c
c     Restore estimates for the normalized derivatives.
c
c-----------------------------------------------------------------------
c
         call dcopy ( n, Work(1,j), 1, Dense(1,j+2), 1 )
c
      end do
c
c-----------------------------------------------------------------------
c
c     Compute coefficients of the Hermite interpolation polynomial.
c
c-----------------------------------------------------------------------
c
      do j = 1, k + 1
         call daxpy ( n, -one, Dense(1,j), 1, Dense(1,j+1), 1 )
      end do
c
c-----------------------------------------------------------------------
c
      return
      end
      subroutine Eval_Herm ( n, Max_Size, Dense, k, tFac, y_Interp )
c
      implicit none
c
c-----------------------------------------------------------------------
c
c     This subroutine evaluates a Hermite interpolation polynomial.
c
c-----------------------------------------------------------------------
c
c     Copyright (C) 2000, Konrad-Zuse-Zentrum fuer Informationstechnik
c     Berlin (ZIB)
c     ALL RIGHTS RESERVED
c
c     Written by:
c
c     R. Ehrig, U. Nowak,
c     Konrad-Zuse-Zentrum fuer Informationstechnik Berlin (ZIB)
c     Takustrasse 7
c     D-14195 Berlin-Dahlem
c
c     phone : +49-30-84185-0
c     fax   : +49-30-84185-125
c     e-mail: ehrig@zib.de, nowak@zib.de
c     URL   : http://www.zib.de
c
c-----------------------------------------------------------------------
c
c     Arguments:
c
c-----------------------------------------------------------------------
c
c     n          Integer variable, must be set by caller on input, not
c                modified. Size of the system.
c
c     Max_Size   Integer variable, must be set by caller on input, not
c                modified. Leading  dimension of  the  array  Dense as
c                declared in LIMEX.
c
c     Dense      Real array  of size at least  Max_Size * k+1, must be
c                set by caller on input. The coefficients of a Hermite
c                interpolation poynomial.
c
c     k          Integer variable, must be set by caller on input. The
c                currrent order of the extrapolation scheme.
c
c     tFac       Real  variable, must  be set  by caller on input. The
c                ratio between t - t0 and the actual stepsize, if t is
c                the current dense output point within this step.
c
c     y_Interp   Real  array of  size n, need not be  set by caller on
c                input. On  exit the estimated  solution values at the
c                dense output  point, i.e. the values  of the  Hermite
c                interpolation polynomials.
c
c-----------------------------------------------------------------------
c
c     Subroutines and functions called
c
c-----------------------------------------------------------------------
c
c     From the BLAS library:
c     ----------------------
c
c     dcopy      Copies real vectors
c
c-----------------------------------------------------------------------
c
c     Declaration of the arguments
c
c-----------------------------------------------------------------------
c
      integer            n, Max_Size, k
c
      double precision   Dense ( Max_Size, * ), tFac, y_Interp ( * )
c
c-----------------------------------------------------------------------
c
c     Local scalars
c
c-----------------------------------------------------------------------
c
      integer            i, j
c
      double precision   tmp
c
c-----------------------------------------------------------------------
c
c     Internal parameters
c
c-----------------------------------------------------------------------
c
      double precision   one
c
      parameter ( one = 1.0d0 )
c
c-----------------------------------------------------------------------
c
c     Evaluate  the Hermite  interpolation polynomials  by the  Horner
c     scheme.
c
c-----------------------------------------------------------------------
c
      tmp = tFac - one
c
      call dcopy ( n, Dense(1,k+2), 1, y_Interp, 1 )
c
      do j = 1, k
         do i = 1, n
            y_Interp(i) = Dense(i,k+2-j) + tmp * y_Interp(i)
         end do
      end do
c
      do i = 1, n
         y_Interp(i) = Dense(i,1) + tFac * y_Interp(i)
      end do
c
c-----------------------------------------------------------------------
c
      return
      end