]>
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 DSPVD2(F,KX,KY,MX,MY,TX,TY,C,NDIMC,NERR) | |
11 | ||
12 | #include "gen/imp64.inc" | |
13 | DIMENSION TX(*),TY(*),C(NDIMC,*) | |
14 | CHARACTER NAME*(*) | |
15 | CHARACTER*80 ERRTXT | |
16 | PARAMETER (NAME = 'DSPVD2') | |
17 | ||
18 | ************************************************************************ | |
19 | * NORBAS, VERSION: 15.03.1993 | |
20 | ************************************************************************ | |
21 | * | |
22 | * DSPVD2 COMPUTES THE COEFFICIENTS | |
23 | * C(I,J) (I=1,...,MX-KX-1 , J=1,...,MY-KY-1) | |
24 | * OF A TWO-DIMENSIONAL POLYNOMIAL VARIATION DIMINISHING SPLINE | |
25 | * APPROXIMATION S(X,Y) IN REPRESENTATION OF NORMALIZED TWO-DIMENSIONAL | |
26 | * B-SPLINES B(I,J)(X,Y) | |
27 | * | |
28 | * S(X,Y) = SUMME(I=1,...,MX-KX-1) | |
29 | * SUMME(J=1,...,MY-KY-1) C(I,J) * B(I,J)(X,Y) . | |
30 | * | |
31 | * THE TWO-DIMENSIONAL B-SPLINES B(I,J)(X,Y) ARE THE PRODUCT OF TWO | |
32 | * ONE-DIMENSIONAL B-SPLINES BX , BY | |
33 | * B(I,J)(X,Y) = BX(I,KX)(X) * BY(J,KY)(Y) | |
34 | * OF DEGREE KX AND KY ( 0 <= KX , KY <= 25 ) WITH INDICES I , J | |
35 | * ( 1 <= I <= MX-KX-1 , 1 <= J <= MY-KY-1 ) OVER TWO SETS OF SPLINE- | |
36 | * KNOTS | |
37 | * TX(1),TX(2),...,TX(MX) ( MX >= 2*KX+2 ) | |
38 | * TY(1),TY(2),...,TY(MY) ( MY >= 2*KY+2 ) , | |
39 | * RESPECTIVELY. | |
40 | * FOR FURTHER DETAILS TO THE ONE-DIMENSIONAL NORMALIZED B-SPLINES SEE | |
41 | * THE COMMENTS TO DSPNB1. | |
42 | * | |
43 | * PARAMETERS: | |
44 | * | |
45 | * F (DOUBLE PRECISION) USER SUPPLIED FUNCTION Z=F(X,Y) FOR WHICH | |
46 | * THE CORRESPONDING SPLINE APPROXIMATION HAS TO BE COMPUTED. | |
47 | * KX (INTEGER) DEGREE OF ONE-DIMENSIONAL B-SPLINES IN X-DIRECTION | |
48 | * OVER THE SET OF KNOTS TX. | |
49 | * KY (INTEGER) DEGREE OF ONE-DIMENSIONAL B-SPLINES IN Y-DIRECTION | |
50 | * OVER THE SET OF KNOTS TY. | |
51 | * MX (INTEGER) NUMBER OF KNOTS FOR THE B-SPLINES IN X-DIRECTION. | |
52 | * MY (INTEGER) NUMBER OF KNOTS FOR THE B-SPLINES IN Y-DIRECTION. | |
53 | * NDIMC (INTEGER) DECLARED FIRST DIMENSION OF ARRAY C IN THE | |
54 | * CALLING PROGRAM, WITH NDIMC >= MX-KX-1 . | |
55 | * TX (DOUBLE PRECISION) ARRAY OF AT LEAST ORDER MX CONTAINING THE | |
56 | * KNOTS IN X-DIRECTION, ON ENTRY. | |
57 | * TY (DOUBLE PRECISION) ARRAY OF AT LEAST ORDER MY CONTAINING THE | |
58 | * KNOTS IN Y-DIRECTION, ON ENTRY. | |
59 | * C (DOUBLE PRECISION) ARRAY OF ORDER (NDIMC, >= MY-KY-1). | |
60 | * ON EXIT C(I,J) CONTAINS THE (I,J)-TH COEFFICIENT OF THE | |
61 | * TWO-DIMENSIONAL B-SPLINE REPRESENTATION OF S(X,Y) . | |
62 | * NERR (INTEGER) ERROR INDICATOR. ON EXIT: | |
63 | * = 0: NO ERROR DETECTED | |
64 | * = 1: AT LEAST ONE OF THE CONSTANTS KX, KY, MX, MY IS ILLEGAL | |
65 | * | |
66 | * ERROR MESSAGES: | |
67 | * | |
68 | * IF ONE OF THE FOLLOWING RELATION IS SATISFIED BY THE CHOSEN INPUT- | |
69 | * PARAMETERS THE PROGRAM RETURNS, AND AN ERROR MESSAGE IS PRINTED: | |
70 | * KX < 1 OR KX > 25 OR KY < 1 OR KY > 25 OR | |
71 | * MX < 2*KX+2 OR MY < 2*KY+2 . | |
72 | * | |
73 | ************************************************************************ | |
74 | ||
75 | PARAMETER (Z1 = 1) | |
76 | ||
77 | NERR=1 | |
78 | IF(KX .LT. 1 .OR. KX .GT. 25) THEN | |
79 | WRITE(ERRTXT,101) 'KX',KX | |
80 | CALL MTLPRT(NAME,'E210.1',ERRTXT) | |
81 | ELSEIF(KY .LT. 1 .OR. KY .GT. 25) THEN | |
82 | WRITE(ERRTXT,101) 'KY',KY | |
83 | CALL MTLPRT(NAME,'E210.1',ERRTXT) | |
84 | ELSEIF(MX .LT. 2*KX+2) THEN | |
85 | WRITE(ERRTXT,101) 'MX',MX | |
86 | CALL MTLPRT(NAME,'E210.2',ERRTXT) | |
87 | ELSEIF(MY .LT. 2*KY+2) THEN | |
88 | WRITE(ERRTXT,101) 'MY',MY | |
89 | CALL MTLPRT(NAME,'E210.2',ERRTXT) | |
90 | ELSE | |
91 | ||
92 | NERR=0 | |
93 | RX=Z1/KX | |
94 | RY=Z1/KY | |
95 | DO 10 I = 1,MX-KX-1 | |
96 | XI=RX*DVSUM(KX,TX(I+1),TX(I+2)) | |
97 | DO 10 J = 1,MY-KY-1 | |
98 | 10 C(I,J)=F(XI,RY*DVSUM(KY,TY(J+1),TY(J+2))) | |
99 | ENDIF | |
100 | ||
101 | RETURN | |
102 | ||
103 | 101 FORMAT(1X,A5,' =',I6,' NOT IN RANGE') | |
104 | END | |
105 | ||
106 | ||
107 |