intrv.f 3.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117
  1. *DECK INTRV
  2. SUBROUTINE INTRV (XT, LXT, X, ILO, ILEFT, MFLAG)
  3. C***BEGIN PROLOGUE INTRV
  4. C***PURPOSE Compute the largest integer ILEFT in 1 .LE. ILEFT .LE. LXT
  5. C such that XT(ILEFT) .LE. X where XT(*) is a subdivision
  6. C of the X interval.
  7. C***LIBRARY SLATEC
  8. C***CATEGORY E3, K6
  9. C***TYPE SINGLE PRECISION (INTRV-S, DINTRV-D)
  10. C***KEYWORDS B-SPLINE, DATA FITTING, INTERPOLATION, SPLINES
  11. C***AUTHOR Amos, D. E., (SNLA)
  12. C***DESCRIPTION
  13. C
  14. C Written by Carl de Boor and modified by D. E. Amos
  15. C
  16. C Abstract
  17. C INTRV is the INTERV routine of the reference.
  18. C
  19. C INTRV computes the largest integer ILEFT in 1 .LE. ILEFT .LE.
  20. C LXT such that XT(ILEFT) .LE. X where XT(*) is a subdivision of
  21. C the X interval. Precisely,
  22. C
  23. C X .LT. XT(1) 1 -1
  24. C if XT(I) .LE. X .LT. XT(I+1) then ILEFT=I , MFLAG=0
  25. C XT(LXT) .LE. X LXT 1,
  26. C
  27. C That is, when multiplicities are present in the break point
  28. C to the left of X, the largest index is taken for ILEFT.
  29. C
  30. C Description of Arguments
  31. C Input
  32. C XT - XT is a knot or break point vector of length LXT
  33. C LXT - length of the XT vector
  34. C X - argument
  35. C ILO - an initialization parameter which must be set
  36. C to 1 the first time the spline array XT is
  37. C processed by INTRV.
  38. C
  39. C Output
  40. C ILO - ILO contains information for efficient process-
  41. C ing after the initial call, and ILO must not be
  42. C changed by the user. Distinct splines require
  43. C distinct ILO parameters.
  44. C ILEFT - largest integer satisfying XT(ILEFT) .LE. X
  45. C MFLAG - signals when X lies out of bounds
  46. C
  47. C Error Conditions
  48. C None
  49. C
  50. C***REFERENCES Carl de Boor, Package for calculating with B-splines,
  51. C SIAM Journal on Numerical Analysis 14, 3 (June 1977),
  52. C pp. 441-472.
  53. C***ROUTINES CALLED (NONE)
  54. C***REVISION HISTORY (YYMMDD)
  55. C 800901 DATE WRITTEN
  56. C 890831 Modified array declarations. (WRB)
  57. C 890831 REVISION DATE from Version 3.2
  58. C 891214 Prologue converted to Version 4.0 format. (BAB)
  59. C 920501 Reformatted the REFERENCES section. (WRB)
  60. C***END PROLOGUE INTRV
  61. C
  62. INTEGER IHI, ILEFT, ILO, ISTEP, LXT, MFLAG, MIDDLE
  63. REAL X, XT
  64. DIMENSION XT(*)
  65. C***FIRST EXECUTABLE STATEMENT INTRV
  66. IHI = ILO + 1
  67. IF (IHI.LT.LXT) GO TO 10
  68. IF (X.GE.XT(LXT)) GO TO 110
  69. IF (LXT.LE.1) GO TO 90
  70. ILO = LXT - 1
  71. IHI = LXT
  72. C
  73. 10 IF (X.GE.XT(IHI)) GO TO 40
  74. IF (X.GE.XT(ILO)) GO TO 100
  75. C
  76. C *** NOW X .LT. XT(IHI) . FIND LOWER BOUND
  77. ISTEP = 1
  78. 20 IHI = ILO
  79. ILO = IHI - ISTEP
  80. IF (ILO.LE.1) GO TO 30
  81. IF (X.GE.XT(ILO)) GO TO 70
  82. ISTEP = ISTEP*2
  83. GO TO 20
  84. 30 ILO = 1
  85. IF (X.LT.XT(1)) GO TO 90
  86. GO TO 70
  87. C *** NOW X .GE. XT(ILO) . FIND UPPER BOUND
  88. 40 ISTEP = 1
  89. 50 ILO = IHI
  90. IHI = ILO + ISTEP
  91. IF (IHI.GE.LXT) GO TO 60
  92. IF (X.LT.XT(IHI)) GO TO 70
  93. ISTEP = ISTEP*2
  94. GO TO 50
  95. 60 IF (X.GE.XT(LXT)) GO TO 110
  96. IHI = LXT
  97. C
  98. C *** NOW XT(ILO) .LE. X .LT. XT(IHI) . NARROW THE INTERVAL
  99. 70 MIDDLE = (ILO+IHI)/2
  100. IF (MIDDLE.EQ.ILO) GO TO 100
  101. C NOTE. IT IS ASSUMED THAT MIDDLE = ILO IN CASE IHI = ILO+1
  102. IF (X.LT.XT(MIDDLE)) GO TO 80
  103. ILO = MIDDLE
  104. GO TO 70
  105. 80 IHI = MIDDLE
  106. GO TO 70
  107. C *** SET OUTPUT AND RETURN
  108. 90 MFLAG = -1
  109. ILEFT = 1
  110. RETURN
  111. 100 MFLAG = 0
  112. ILEFT = ILO
  113. RETURN
  114. 110 MFLAG = 1
  115. ILEFT = LXT
  116. RETURN
  117. END