]>
Commit | Line | Data |
---|---|---|
fe4da5cc | 1 | * |
2 | * $Id$ | |
3 | * | |
4 | * $Log$ | |
5 | * Revision 1.1.1.1 1996/04/01 15:02:25 mclareni | |
6 | * Mathlib gen | |
7 | * | |
8 | * | |
9 | #include "gen/pilot.h" | |
10 | SUBROUTINE DSPAP1(K,M,N,XI,YI,KNOT,T,C,W,NW,NERR) | |
11 | ||
12 | #include "gen/imp64.inc" | |
13 | DIMENSION XI(*),YI(*),T(*),W(*),C(*) | |
14 | CHARACTER NAME*(*) | |
15 | CHARACTER*80 ERRTXT | |
16 | PARAMETER (NAME = 'DSPAP1') | |
17 | ||
18 | ************************************************************************ | |
19 | * NORBAS, VERSION: 15.03.1993 | |
20 | ************************************************************************ | |
21 | * | |
22 | * DSPAP1 COMPUTES THE COEFFICIENTS C(1),...,C(NC) OF A POLYNOMIAL | |
23 | * APPROXIMATION SPLINE S(X) IN B-SPLINE REPRESENTATION | |
24 | * | |
25 | * S(X) = SUMME(I=1,...,NC) C(I) * B(I,K)(X) , NC=M-K-1 | |
26 | * | |
27 | * TO A USER SUPPLIED DATA SET | |
28 | * | |
29 | * (XI(J),YI(J)) , J = 1,2,...,N , N >= M-K-1 >= K+1 | |
30 | * | |
31 | * OF A FUNCTION Y=F(X) , I.E. | |
32 | * | |
33 | * S(XI(J)) = YI(J) , J = 1,2,...,N . | |
34 | * | |
35 | * THE FUNCTIONS B(I,K)(X) ARE NORMALIZED B-SPLINES OF DEGREE K | |
36 | * (0<= K <= 25) WITH INDEX I (1 <= I <= N) OVER A SET OF SPLINE-KNOTS | |
37 | * T(1),T(2), ... ,T(M) ( M <= N+K+1 ) | |
38 | * (KNOTS IN ASCENDING ORDER, WITH MULTIPLICITIES NOT GREATER | |
39 | * THAN K+1). | |
40 | * FOR FURTHER DETAILS TO THE ONE-DIMENSIONAL NORMALIZED B-SPLINES SEE | |
41 | * THE COMMENTS TO DSPNB1. | |
42 | * | |
43 | * PARAMETERS: | |
44 | * | |
45 | * N (INTEGER) NUMBER OF APPROXIMATION POINTS . | |
46 | * M (INTEGER) NUMBER OF KNOTS. | |
47 | * K (INTEGER) DEGREE OF B-SPLINES. | |
48 | * XI (DOUBLE PRECISION) ARRAY OF AT LEAST ORDER N . | |
49 | * XI MUST CONTAIN THE APPROXIMATION POINTS IN ASCENDING ORDER, | |
50 | * ON ENTRY. | |
51 | * YI (DOUBLE PRECISION) ARRAY OF AT LEAST ORDER N CONTAINING | |
52 | * THE FUNCTION VALUES YI(J), J=1,...,N, ON ENTRY. | |
53 | * KNOT (INTEGER) PARAMETER FOR STEERING THE CHOICE OF KNOTS. | |
54 | * ON ENTRY: | |
55 | * = 1 : KNOTS ARE COMPUTED BY DSPAP1 IN THE FOLLOWING WAY: | |
56 | * T(J) = XI(1) , J = 1,...,K+1 | |
57 | * T(J) = XI(1)+(J-K-1)*(XI(N)-XI(1))/(NC-K) , | |
58 | * J = K+2,...,NC | |
59 | * T(NC+J) = XI(N) , J = 1,...,K+1 | |
60 | * = 2 : KNOTS ARE COMPUTED BY DSPAP1 IN THE FOLLOWING WAY: | |
61 | * T(J) = XI(1) , J = 1,...,K+1 | |
62 | * T(J) = (XI(J-K-1)+XI(J))/2 , J = K+2,...,N | |
63 | * T(N+J) = XI(N) , J = 1,...,K+1 | |
64 | * OTHERWISE KNOTS ARE USER SUPPLIED. RECOMMENDED CHOICE : | |
65 | * T(1) <= ... <= T(K+1) <= XI(1) | |
66 | * XI(1) < T(K+2) < ... < T(NC) < XI(N) | |
67 | * XI(N) <= T(NC+1) <= ... <= T(NC+K+1) | |
68 | * T (DOUBLE PRECISION) ARRAY OF AT LEAST ORDER M . | |
69 | * IF THE INPUT VALUE OF THE PARAMETER KNOT IS 1 OR 2 THE | |
70 | * KNOTS ARE COMPUTED BY DSPAP1 AND THEY ARE GIVEN IN THE | |
71 | * ARRAY T, ON EXIT. | |
72 | * IN THE OTHER CASES THE ARRAY T MUST CONTAIN THE USER | |
73 | * SUPPLIED KNOTS, ON ENTRY. | |
74 | * W (DOUBLE PECISION) WORKING ARRAY OF AT LEAST ORDER NW. | |
75 | * NW (INTEGER) ORDER OF WORKING ARRAY W . | |
76 | * NW >= N*(NC+5)+NC*(NC+1) , WITH NC=M-K-1 . | |
77 | * FOR GOOD PERFORMANCE, NW SHOULD GENERALLY BE LARGER. | |
78 | * C (DOUBLE PRECISION) ARRAY OF AT LEAST ORDER N . ON EXIT | |
79 | * C(1),...,C(NC) CONTAIN THE COEFFICIENTS OF THE B-SPLINE | |
80 | * REPRESENTATION OF S(X). | |
81 | * NERR (INTEGER) ERROR INDICATOR. ON EXIT: | |
82 | * = 0: NO ERROR DETECTED | |
83 | * = 1: AT LEAST ONE OF THE CONSTANTS K , M , N IS ILLEGAL | |
84 | * | |
85 | * ERROR MESSAGES: | |
86 | * | |
87 | * IF ONE OF THE FOLLOWING RELATION IS SATISFIED BY THE CHOSEN INPUT- | |
88 | * PARAMETERS THE PROGRAM RETURNS, AND AN ERROR MESSAGE IS PRINTED: | |
89 | * K < 0 OR K > 25 OR | |
90 | * M < 2*K+2 OR | |
91 | * N < M-K-1 . | |
92 | * | |
93 | ************************************************************************ | |
94 | ||
95 | NERR=1 | |
96 | IF(K .LT. 0 .OR. K .GT. 25) THEN | |
97 | WRITE(ERRTXT,101) 'K',K | |
98 | CALL MTLPRT(NAME,'E210.1',ERRTXT) | |
99 | ELSEIF(M .LT. 2*K+2) THEN | |
100 | WRITE(ERRTXT,101) 'M',M | |
101 | CALL MTLPRT(NAME,'E210.2',ERRTXT) | |
102 | ELSEIF(N .LT. M-K-1) THEN | |
103 | WRITE(ERRTXT,101) 'N',N | |
104 | CALL MTLPRT(NAME,'E210.4',ERRTXT) | |
105 | ELSE | |
106 | ||
107 | NC=M-K-1 | |
108 | M1=1 | |
109 | M2=M1+N*NC | |
110 | M3=M2+NC | |
111 | M4=M3+NC*NC | |
112 | LW=NW-M4+1 | |
113 | ||
114 | CALL SPLAS1(N,NC,M,K,XI,YI,KNOT,T,W,W(M2),W(M3),W(M4),LW,C,NERR) | |
115 | ||
116 | ENDIF | |
117 | RETURN | |
118 | ||
119 | 101 FORMAT(1X,A5,' =',I6,' NOT IN RANGE') | |
120 | END | |
121 | ||
122 | ||
123 |