interp_linear Subroutine

subroutine interp_linear ( dim_num, data_num, t_data, p_data, interp_num, &
     t_interp, p_interp )

  !*****************************************************************************80
  !
  !! INTERP_LINEAR applies piecewise linear interpolation to data.
  !
  !  Discussion:
  !
  !    From a space of DIM_NUM dimensions, we are given a sequence of
  !    DATA_NUM points, which are presumed to be successive samples
  !    from a curve of points P.
  !
  !    We are also given a parameterization of this data, that is,
  !    an associated sequence of DATA_NUM values of a variable T.
  !    The values of T are assumed to be strictly increasing.
  !
  !    Thus, we have a sequence of values P(T), where T is a scalar,
  !    and each value of P is of dimension DIM_NUM.
  !
  !    We are then given INTERP_NUM values of T, for which values P
  !    are to be produced, by linear interpolation of the data we are given.
  !
  !    Note that the user may request extrapolation.  This occurs whenever
  !    a T_INTERP value is less than the minimum T_DATA or greater than the
  !    maximum T_DATA.  In that case, linear extrapolation is used.  
  !
  !  Licensing:
  !
  !    This code is distributed under the GNU LGPL license. 
  !
  !  Modified:
  !
  !    03 December 2007
  !
  !  Author:
  !
  !    John Burkardt
  !
  !  Parameters:
  !
  !    Input, integer ( kind = 4 ) DIM_NUM, the spatial dimension.
  !
  !    Input, integer ( kind = 4 ) DATA_NUM, the number of data points.
  !
  !    Input, real ( kind = 8 ) T_DATA(DATA_NUM), the value of the
  !    independent variable at the sample points.  The values of T_DATA
  !    must be strictly increasing.
  !
  !    Input, real ( kind = 8 ) P_DATA(DIM_NUM,DATA_NUM), the value of the
  !    dependent variables at the sample points.
  !
  !    Input, integer ( kind = 4 ) INTERP_NUM, the number of points
  !    at which interpolation is to be done.
  !
  !    Input, real ( kind = 8 ) T_INTERP(INTERP_NUM), the value of the
  !    independent variable at the interpolation points.
  !
  !    Output, real ( kind = 8 ) P_INTERP(DIM_NUM,DATA_NUM), the interpolated
  !    values of the dependent variables at the interpolation points.
  !
  implicit none

  integer ( kind = 4 ) data_num
  integer ( kind = 4 ) dim_num
  integer ( kind = 4 ) interp_num

  integer ( kind = 4 ) interp
  integer ( kind = 4 ) left
  real    ( kind = 8 ) p_data(dim_num,data_num)
  real    ( kind = 8 ) p_interp(dim_num,interp_num)
  logical              r8vec_ascends_strictly
  integer ( kind = 4 ) right
  real    ( kind = 8 ) t
  real    ( kind = 8 ) t_data(data_num)
  real    ( kind = 8 ) t_interp(interp_num)

  if ( .not. r8vec_ascends_strictly ( data_num, t_data ) ) then
     write ( *, '(a)' ) ' '
     write ( *, '(a)' ) 'INTERP_LINEAR - Fatal error!'
     write ( *, '(a)' ) '  Independent variable array T_DATA is not strictly increasing.'
     stop
  end if

  do interp = 1, interp_num

     t = t_interp(interp)
     !
     !  Find the interval [ TDATA(LEFT), TDATA(RIGHT) ] that contains, or is
     !  nearest to, TVAL.
     !
     call r8vec_bracket ( data_num, t_data, t, left, right )

     p_interp(1:dim_num,interp) = &
          ( ( t_data(right) - t                ) * p_data(1:dim_num,left)   &
          + (                 t - t_data(left) ) * p_data(1:dim_num,right) ) &
          / ( t_data(right)     - t_data(left) )

  end do

  return
end subroutine interp_linear