1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
|
SUBROUTINE CONINT (NDP,XD,YD,ZD,NCP,IPC,PD)
C
C +-----------------------------------------------------------------+
C | |
C | Copyright (C) 1986 by UCAR |
C | University Corporation for Atmospheric Research |
C | All Rights Reserved |
C | |
C | NCARGRAPHICS Version 1.00 |
C | |
C +-----------------------------------------------------------------+
C
C
C
C THIS SUBROUTINE ESTIMATES PARTIAL DERIVATIVES OF THE FIRST AND
C SECOND ORDER AT THE DATA POINTS.
C THE INPUT PARAMETERS ARE
C
C NDP = NUMBER OF DATA POINTS,
C XD,YD,ZD = ARRAYS CONTAINING THE X, Y, AND Z COORDI-
C NATES OF DATA POINTS,
C NCP = NUMBER OF DATA POINTS TO BE USED FOR ESTIMATION
C OF PARTIAL DERIVATIVES AT EACH DATA POINT,
C IPC = INTEGER ARRAY CONTAINING THE POINT NUMBERS OF
C NCP DATA POINTS CLOSEST TO EACH OF THE NDP DATA
C POINT.
C THE OUTPUT PARAMETER IS
C
C PD = ARRAY OF DIMENSION 5*NDP, WHERE THE ESTIMATED
C
C ZX, ZY, ZXX, ZXY, AND ZYY VALUES AT THE DATA
C POINTS ARE TO BE STORED.
C DECLARATION STATEMENTS
C
C
DIMENSION XD(NDP) ,YD(NDP) ,ZD(NDP) ,IPC(1) ,
1 PD(1)
REAL NMX ,NMY ,NMZ ,NMXX ,
1 NMXY ,NMYX ,NMYY
C
SAVE
C
C PRELIMINARY PROCESSING
C
C
NCPM1 = NCP-1
C
C ESTIMATION OF ZX AND ZY
C
C
DO 130 IP0=1,NDP
X0 = XD(IP0)
Y0 = YD(IP0)
Z0 = ZD(IP0)
NMX = 0.0
NMY = 0.0
NMZ = 0.0
JIPC0 = NCP*(IP0-1)
DO 120 IC1=1,NCPM1
JIPC = JIPC0+IC1
IPI = IPC(JIPC)
DX1 = XD(IPI)-X0
DY1 = YD(IPI)-Y0
DZ1 = ZD(IPI)-Z0
IC2MN = IC1+1
DO 110 IC2=IC2MN,NCP
JIPC = JIPC0+IC2
IPI = IPC(JIPC)
DX2 = XD(IPI)-X0
DY2 = YD(IPI)-Y0
DNMZ = DX1*DY2-DY1*DX2
IF (DNMZ .EQ. 0.0) GO TO 110
DZ2 = ZD(IPI)-Z0
DNMX = DY1*DZ2-DZ1*DY2
DNMY = DZ1*DX2-DX1*DZ2
IF (DNMZ .GE. 0.0) GO TO 100
DNMX = -DNMX
DNMY = -DNMY
DNMZ = -DNMZ
100 NMX = NMX+DNMX
NMY = NMY+DNMY
NMZ = NMZ+DNMZ
110 CONTINUE
120 CONTINUE
JPD0 = 5*IP0
PD(JPD0-4) = -NMX/NMZ
PD(JPD0-3) = -NMY/NMZ
130 CONTINUE
C
C ESTIMATION OF ZXX, ZXY, AND ZYY
C
C
DO 170 IP0=1,NDP
JPD0 = JPD0+5
X0 = XD(IP0)
JPD0 = 5*IP0
Y0 = YD(IP0)
ZX0 = PD(JPD0-4)
ZY0 = PD(JPD0-3)
NMXX = 0.0
NMXY = 0.0
NMYX = 0.0
NMYY = 0.0
NMZ = 0.0
JIPC0 = NCP*(IP0-1)
DO 160 IC1=1,NCPM1
JIPC = JIPC0+IC1
IPI = IPC(JIPC)
DX1 = XD(IPI)-X0
DY1 = YD(IPI)-Y0
JPD = 5*IPI
DZX1 = PD(JPD-4)-ZX0
DZY1 = PD(JPD-3)-ZY0
IC2MN = IC1+1
DO 150 IC2=IC2MN,NCP
JIPC = JIPC0+IC2
IPI = IPC(JIPC)
DX2 = XD(IPI)-X0
DY2 = YD(IPI)-Y0
DNMZ = DX1*DY2-DY1*DX2
IF (DNMZ .EQ. 0.0) GO TO 150
JPD = 5*IPI
DZX2 = PD(JPD-4)-ZX0
DZY2 = PD(JPD-3)-ZY0
DNMXX = DY1*DZX2-DZX1*DY2
DNMXY = DZX1*DX2-DX1*DZX2
DNMYX = DY1*DZY2-DZY1*DY2
DNMYY = DZY1*DX2-DX1*DZY2
IF (DNMZ .GE. 0.0) GO TO 140
DNMXX = -DNMXX
DNMXY = -DNMXY
DNMYX = -DNMYX
DNMYY = -DNMYY
DNMZ = -DNMZ
140 NMXX = NMXX+DNMXX
NMXY = NMXY+DNMXY
NMYX = NMYX+DNMYX
NMYY = NMYY+DNMYY
NMZ = NMZ+DNMZ
150 CONTINUE
160 CONTINUE
PD(JPD0-2) = -NMXX/NMZ
PD(JPD0-1) = -(NMXY+NMYX)/(2.0*NMZ)
PD(JPD0) = -NMYY/NMZ
170 CONTINUE
RETURN
END
|
