Search Site | Contact Details | FAQ

ADAS Subroutine xxnbaf

      SUBROUTINE XXNBAF (M, NCAP7 , X, Y, W, LAMDA , B, A, DIAG, C ,
     &                   SS , IFAIL )
C
C
C-----------------------------------------------------------------------
C
C  PURPOSE: Determines a least-square cubic spline approximation s(x)
C           to the set of data points (x_r, y_r) with weights w_r,
C	    for r=1,2,...,m.
C
C           The value of NCAP7 = ncap+7, where ncap is the number of
C	    intervals of the spline (number of interior knots + 1),
C           and the values of the knots LAMDA(5), LAMDA(6), ...,
C	    LAMDA(NCAP7-4), interior to the data interval, are
C	    prescribed by the user.
C
C	    s has the property that it minimises ss, the sum of the
C	    squares of the weighted residuals eps(r)
C
C		      eps(r) = w(r) * ( s(x(r)) - y(r) ).
C
C	    The procedure produces the minimising value of ss and
C	    the coefficients c(1), c(2),...,c(q), where q=ncap+3=NCAP7-4,
C	    in the B-spline representation
C
C		s(x) = c(1)*N1(x) + c(2)*N2(x) + ... + c(q)*Nq(x) .
C
C	    Here Ni(x) (i=1,2,...,q) denotes the normalised B-spline
C	    of degree 3 defined upon the knots lamda(i-4), lamda(i-3),
C	    lamda(i-2), lamda(i-1), and lamda(i).
C
C  CALLING PROGRAM: VARIOUS
C
C  SUBROUTINE:
C
C  INPUT: (I*4) M            = The number of data points.
C                              CONSTRAINT: M >= MDIST >= 4, where
C                              MDIST is the number of distinct x
C                              values in the data.
C
C  INPUT: (I*4) NCAP7        = NBAR+7, where NBAR is the number of
C                              intervals of the spline (number
C                              of interior knots +1, i.e. the knots 
C                              strictly in the range X(1) to X(M)) 
C                              over which the spline is defined.
C                              CONSTRAINT: 8<= NCAP7 <= MDIST+4,
C                              where MDIST is the number of distinct
C                              x values in the data. 
C
C  INPUT:  (R*8) X()         = The values x_r of the independent variable 
C                              (abscissa), for r=1,2,...,m.
C                              CONSTRAINT: 
C                                 X(1)<=X(2)<=...<=X(M)          
C
C  INPUT:  (R*8) Y()         = The values y_r of the dependent variable 
C                              (ordinate), for r=1,2,...,m.
C
C  INPUT:  (R*8) W()         = The values w_r of the weights,           
C                              for r=1,2,...,m.
C
C  INPUT:  (R*8) LAMDA()     = LAMDA(i) must be set to the (i-4)th
C                              (interior) knot, i=5,6,...,nbar+3.
C                              CONSTRAINT:                          
C                                 X(1) < LAMDA(5) <= LAMDA(6) ... <=
C                                       ... <=LAMDA(NCAP7-4) < X(M) .
C
C  INPUT:  (I*8) IFAIL       = 0 : stop if any error
C                            = 1 : continue if non-fatal error.
C
C  OUTPUT: (R*8) LAMDA()     = Input values are unchanged, and
C                              LAMDA(i), for i=1,2,3,4,NCAP7-3,
C                              NCAP7-2,NCAP7-1,NCAP7 contains the
C                              additional exterior knots introduced by
C                              the routine.
C
C
C  OUTPUT: (R*8) C()         = The coefficients of the B-spline N_i(x),
C                              for i=1,2,...,nbar+3. The remaining 
C                              elements (from NBAR+4 to NBAR+7) are not
C                              used.
C
C  OUTPUT: (R*8) SS          = The residual sum of sqaures
C
C  OUTPUT: (I*4) IFAIL       = 0 : no error detected
C                            = 1 : the knots fail to satisfy the condition
C                                   X(1) < LAMDA(5) <= LAMDA(6) <=...
C                                     <= LAMDA(NCAP7-4) < X(M)
C                            = 2 : The weights are not strictly positive
C                            = 3 : The values of X(R), R=1,M are not in
C                                  non-decreasing order.
C                            = 4 : NCAP7 < 8 (so that the number of
C                                  interior knots is negative) or
C                                  NCAP7 > MDIST + 4, where MDIST is the 
C                                  number of distinct x values in the data
C                                  (so that there cannot be unique solution).
C                            = 5 : The conditions specified by Schoenberg
C                                  and Whitney fail to hold for at least
C                                  one subset of the distinct data abscissae.
C                                  That is, there is no subset of NCAP7-4
C                                  strictly increasing values,
C                                  X(R(1)), X(R(2)),..., X(R(NCAP7-4)),
C                                  among theabscissae such that
C
C                                  X(R(1)) < LAMDA(1) < X(R(5)) 
C                                  X(R(2)) < LAMDA(2) < X(R(6)) 
C                                  ...
C                                  X(R(NCAP7-8))<LAMDA(NCAP7-8)<X(R(NCAP7-4)).
C
C                                  This means that there is no unique 
C                                  solution: there are regions containing
C                                  too many knots compared with the
C                                  number of data points.
C
C          (R*8)  B()        = Set of distinct data abscissae
C    
C          (R*8)  WORK2()    = WORKSPACE
C
C          (I*4)  J          = GENERAL INDEX
C          (I*4)  I          = GENERAL INDEX
C          (I*4)  R          = GENERAL INDEX
C          (I*4)  II         = GENERAL INDEX
C          (R*8)  BI         = GENERAL REAL
C          (R*8)  XI         = GENERAL REAL
C
C  ROUTINES:  NONE
C
C  AUTHORS: Alessandro C. Lanzafame, University of Strathclyde
C
C  REFERENCE: Cox, M.G. and Hayes, J.G. "Curve fitting: A Guide and
C             Suite of Algorithms for the Non-specialist User."
C             Report NAC26, National Physical Laboratory, Middlessex,
C             1973.
C
C  DATE:    12 January 1995
C
C  VERSIION: 1.0a
C  Alessandro Lanzafame, 12 January 1995.
C  Directly derived from Algol text.
C  (Error in passing woking variables)
C
C  VERSION 1.0b
C  Alessandro Lanzafame, 15 January 1995.
C  DIAG(1:NCAP7-4) absorbed in matrix A(1:NCAP7-4,1). 
C  Matrix A(1:NCAP7-4,2:4) becomes A(1:NCAP7-4,1:4). This is to easy
C  the passing of workspaces.
C  WORK1 identified with B. WORK2 with A.
C  Corrected index error in checking remaining Schoennber-Whitney 
C  conditions.
C  COMPILING BUT NOT WORKING.
C
C  VERSION 1.0c
C  Alessandro Lanzafame, 15 January 1995.
C  Knots shifted as in original Algol program text.
C
C UNIX-IDL PORT:
C
C VERSION: 1.1                          DATE: 22-1-96
C MODIFIED: TIM HAMMOND (TESSELLA SUPPORT SERVICES PLC)
C               - PUT UNDER SCCS CONTROL
C
C VERSION: 1.2                          DATE: 06-07-2004
C MODIFIED: Allan Whiteford
C               - Changed name from dxnbaf to xxnbaf.
C
C VERSION  : 1.3			DATE: 10-04-2007
C MODIFIED : Allan Whiteford
C               - Modified documentation as part of automated
C		  subroutine documentation preparation.
C
C-----------------------------------------------------------------------
C
C-----------------------------------------------------------------------
      INTEGER             IFAIL,       M,           NCAP7
      REAL*8              A(1:NCAP7-4,2:4),         B(M)
      REAL*8              C(NCAP7),    DIAG(1:NCAP7-4)
      REAL*8              LAMDA(-3:NCAP7-4),        SS,          W(M)
      REAL*8              X(M),        Y(M)
© Copyright 1995-2024 The ADAS Project
Comments and questions to: adas-at-adas.ac.uk