]>
Commit | Line | Data |
---|---|---|
fe4da5cc | 1 | * |
2 | * $Id$ | |
3 | * | |
4 | * $Log$ | |
5 | * Revision 1.1.1.1 1996/04/01 15:02:26 mclareni | |
6 | * Mathlib gen | |
7 | * | |
8 | * | |
9 | #include "gen/pilot.h" | |
10 | SUBROUTINE SPLAS1(N,NC,M,K,XI,YI,KNOT,T,A,S,VT,W,LW,C,NERR) | |
11 | ||
12 | #include "gen/imp64.inc" | |
13 | DIMENSION XI(*),YI(*),T(*),A(N,*),S(*),VT(NC,*),W(*),C(*) | |
14 | ||
15 | ************************************************************************ | |
16 | * NORBAS, VERSION: 15.03.1993 | |
17 | ************************************************************************ | |
18 | * | |
19 | * THE SUBROUTINE SPLAS1 IS USED BY DSPAP1 FOR COMPUTING THE | |
20 | * COEFFICIENTS C(1),...,C(NC) OF A POLYNOMIAL APPROXIMATION SPLINE | |
21 | * S(X) IN B-SPLINE REPRESENTATION | |
22 | * | |
23 | ************************************************************************ | |
24 | ||
25 | PARAMETER (Z0 = 0 , Z1 = 1 , Z2 = 2 , Z10 = 10 , HALF = Z1/Z2) | |
26 | * | |
27 | * COMPUTE AN APPROXIMATION EPS0 TO THE RELATIVE MACHINE PRECISION | |
28 | * | |
29 | EPS0=Z1 | |
30 | 10 EPS0=EPS0/Z10 | |
31 | IF (Z1+EPS0 .NE. Z1) GO TO 10 | |
32 | EPS0=Z10*EPS0 | |
33 | * | |
34 | * COMPUTE KNOTS BY MEANS OF GIVEN DATA POINTS (IF KNOT = 1 OR 2) | |
35 | * | |
36 | IF (KNOT .EQ. 1) THEN | |
37 | CALL DSPKN1(K,M,XI(1),XI(N),T,NERR) | |
38 | ELSEIF (KNOT .EQ. 2) THEN | |
39 | DO 20 I=1,K+1 | |
40 | T(I)=XI(1) | |
41 | 20 T(NC+I)=XI(N) | |
42 | DO 30 I=K+2,NC | |
43 | 30 T(I)=HALF*(XI(N*(I-K-2)/NC+1)+XI(N*I/NC)) | |
44 | ENDIF | |
45 | * | |
46 | * COMPUTE MATRIX A AND SOLVE LINEAR LEAST SQUARES PROBLEM USING SVD | |
47 | * | |
48 | DO 40 I=1,N | |
49 | DO 40 J=1,NC | |
50 | 40 A(I,J)=DSPNB1(K,M,J,0,XI(I),T,NERR) | |
51 | CALL DGESVD('O','A',N,NC,A,N,S,W,1,VT,NC,W,LW,INFO) | |
52 | CALL DMMPY(NC,N,A(1,1),A(2,1),A(1,2),YI(1),YI(2),W(1),W(2)) | |
53 | DO 50 J=1,NC | |
54 | IF (S(J) .GT. EPS0*S(1)) THEN | |
55 | W(J)=W(J)/S(J) | |
56 | ELSE | |
57 | W(J)=Z0 | |
58 | ENDIF | |
59 | 50 CONTINUE | |
60 | CALL DMMPY(NC,NC,VT(1,1),VT(2,1),VT(1,2),W(1),W(2),C(1),C(2)) | |
61 | NERR=0 | |
62 | ||
63 | RETURN | |
64 | END | |
65 | ||
66 | ||
67 |