interp1D Subroutine

public subroutine interp1D(tabgros, lb_gros, tabfin, lb_fin, ub_gros, ordre)

Interpolate evenly spaced points, taking into account the borders

Arguments

Type IntentOptional Attributes Name
real(kind=R8), intent(in), dimension(lb_gros:) :: tabgros

tableau grossier à interpoler

integer(kind=I4), intent(in) :: lb_gros

indice inférieur

real(kind=R8), intent(out), dimension(lb_fin :) :: tabfin

tableau résultant, 2 fois plus fin

integer(kind=I4), intent(in) :: lb_fin

indice inférieur de tab_fin

integer(kind=I4), intent(in) :: ub_gros

taille de tabgros

integer(kind=I4), intent(in) :: ordre

ordre de l’interpolation


Calls

proc~~interp1d~~CallsGraph proc~interp1d interp1D proc~interp interp proc~interp1d->proc~interp

Called by

proc~~interp1d~~CalledByGraph proc~interp1d interp1D proc~interp2d interp2D proc~interp2d->proc~interp1d proc~test_interp_pond test_interp_pond proc~test_interp_pond->proc~interp1d proc~test_interp_pond->proc~interp2d program~test_intpl test_intpl program~test_intpl->proc~test_interp_pond

Source Code

subroutine interp1D(tabgros, lb_gros, tabfin, lb_fin, ub_gros, ordre)
!! Interpolate evenly spaced points, taking into account the borders
implicit none
integer(kind=I4), intent(in)                       :: lb_gros  !! *indice inférieur*
integer(kind=I4), intent(in)                       :: lb_fin   !! *indice inférieur de tab_fin*
integer(kind=I4), intent(in)                       :: ub_gros  !! *taille de tabgros*
integer(kind=I4), intent(in)                       :: ordre    !! *ordre de l'interpolation*
real(kind=R8),    intent(in), dimension(lb_gros:)  :: tabgros  !! *tableau grossier à interpoler*
real(kind=R8),   intent(out), dimension(lb_fin :)  :: tabfin   !! *tableau résultant, 2 fois plus fin*

   integer(kind=I4) :: l_inf, l_sup, i, ii
   real(kind=R8)    :: tmp0, dtmp

   real(kind=R8), dimension(       -ordre/2:        ordre  ) :: tab_inf
   real(kind=R8), dimension(ub_gros-ordre  :ub_gros+ordre/2) :: tab_sup

   ! bornes pour déterminer les limites d'utilisation de la fonction interp
   l_inf = ordre/2
   l_sup = ub_gros -l_inf

   ! extension du tableau par prolongement de la dérivée
   tab_inf(0:ordre) = tabgros(0:ordre)
   tmp0 = tab_inf(0)
   dtmp = tab_inf(0)-tab_inf(1)

   do i = 1, l_inf
      tab_inf(-i) = tmp0 + i*dtmp
   enddo

   do ii = 0, l_inf
      i = 2*ii
      tabfin(i  ) = tabgros(ii)
      tabfin(i+1) = interp( tab   = tab_inf,    &  !
                            lb    = -ordre/2,   &  !
                            ind   = ii,         &  !
                            ordre = ordre )        !
   enddo

   ! utilisation d'interp dans les limites normales
   do ii = l_inf+1, l_sup-1
      i = 2*ii
      tabfin(i  ) = tabgros(ii)
      tabfin(i+1) = interp( tab   = tabgros,    &  !
                            lb    = lb_gros,    &  !
                            ind   = ii,         &  !
                            ordre = ordre )        !
   enddo

   ! extension du tableau par prolongement de la dérivée
   tab_sup(ub_gros-ordre:ub_gros) = tabgros(ub_gros-ordre:ub_gros)
   tmp0 = tab_sup(ub_gros)
   dtmp = tab_sup(ub_gros)-tab_sup(ub_gros-1)

   do i = 1, l_inf
      tab_sup(i+ub_gros) = tmp0 + i*dtmp
   enddo

   do ii = l_sup, ub_gros-1
      i = 2*ii
      tabfin(i  ) = tabgros(ii)
      tabfin(i+1) = interp( tab   = tab_sup,          &  !
                            lb    = ub_gros - ordre,  &  !
                            ind   = ii,               &  !
                            ordre = ordre )              !
   enddo
   tabfin(2*ub_gros) = tabgros(ub_gros)

return
endsubroutine interp1D