summaryrefslogtreecommitdiff
path: root/Donjon/src/DETSPL3.f
blob: 732736640becfef929e006eac41452659f7b9971 (plain)
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
*DECK DETSPL3
      SUBROUTINE DETSPL3(XCNTR ,YCNTR ,ZCNTR ,
     >                  NXMAX ,NYMAX ,NZMAX ,
     >                  FXYZ  ,FXY   ,FDUMMY,
     >                  F2X   ,F2Y   ,F2Z   ,
     >                  XINT  ,YINT  ,ZINT  ,
     >                  FP1   ,FP2   ,
     >                  FYINT ,FZINT ,FINTRP,
     >                  N1    ,N2    ,N3    ,ITYPE)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Perform spline interpolation.
*
*Copyright:
* Copyright (C) 2010 Ecole Polytechnique de Montreal.
*
*Author(s): 
* E. Varin
*
*
*Parameters: 
* XCNTR     
* YCNTR     
* ZCNTR     
* NXMAX     
* NYMAX     
* NZMAX     
* FXYZ      
* FXY       
* FDUMMY    
* F2X       
* F2Y       
* F2Z       
* XINT      
* YINT      
* ZINT      
* FP1       
* FP2       
* FYINT     
* FZINT     
* FINTRP    
* N1        
* N2        
* N3        
* ITYPE     
*
*-----------------------------------------------------------------------
*
      REAL*4     XCNTR(NXMAX)  ,YCNTR(NYMAX)     ,ZCNTR(NZMAX),
     >           FXYZ(N1,N2,N3),FXY(NXMAX,NYMAX) ,
     >           FYINT(NYMAX)  ,FZINT(NZMAX)     ,
     >           F2X(NXMAX)    ,F2Y(NYMAX)       ,F2Z(NZMAX),
     >           FDUMMY(NXMAX)
*----
*  INTERPOLATE IN TWO DIMENSIONS AT XINT,YINT FOR EACH Z PLANE
*----
      DO 10 K=1,NZMAX

         DO 20 J=1,NXMAX
            DO 30 I=1,NYMAX

               IF (ITYPE.EQ.1) THEN
                  FXY(J,I) = FXYZ(J,I,K)
               ELSE IF (ITYPE.EQ.2) THEN
                  FXY(J,I) = FXYZ(K,J,I)
               ELSE IF (ITYPE.EQ.3) THEN
                  FXY(J,I) = FXYZ(I,K,J)
               ELSE
                  CALL XABORT('DETSPL3: ERROR IN SPLIN3')
               ENDIF

   30       CONTINUE
  20     CONTINUE

         CALL DETSPL2(XCNTR,YCNTR,NXMAX ,NYMAX ,FXY,
     >               FP1  ,FP2  ,F2X   ,F2Y   ,FDUMMY,
     >               XINT ,YINT ,FYINT ,FXYINT)

         FZINT(K) = FXYINT

   10 CONTINUE
*----
*  CALCULATE SECOND DERIVATIVE ALONG Z AT XINT,YINT
*----
      CALL DETSPLI(ZCNTR,FZINT,NZMAX,FP1,FP2,F2Z)
*----
*  INTERPOLATE ALONG Z FOR XINT,YINT
*----
      CALL DETSPLI2(ZCNTR,FZINT,F2Z,NZMAX,ZINT,FINTRP)

      RETURN
      END