Example 4 -- interpolated parametric curves and derivatives


Please note that to compile the graphics portion of the code below, you must link with a double precision version of NCAR Graphics, or else you need to convert the code back to single precision before you call the graphics routines.

      PROGRAM FTEX04
C
C  Example of KURV1DP, KURV2DP, KURVDDP.
C
      PARAMETER (IDIM=11,IOUT=201,IDTEMP=9*IDIM)
      DOUBLE PRECISION X
      DOUBLE PRECISION Y
      DOUBLE PRECISION TEMP
      DOUBLE PRECISION U
      DOUBLE PRECISION XO
      DOUBLE PRECISION YO
      DOUBLE PRECISION XS
      DOUBLE PRECISION YS
      DOUBLE PRECISION XD
      DOUBLE PRECISION YD
      DOUBLE PRECISION XDD
      DOUBLE PRECISION YDD
      DOUBLE PRECISION XP
      DOUBLE PRECISION YP
      DOUBLE PRECISION S
      DOUBLE PRECISION SIGMA
      DOUBLE PRECISION SLOP1
      DOUBLE PRECISION SLOP2
      DOUBLE PRECISION TINC
      DIMENSION X(IDIM),Y(IDIM),TEMP(IDTEMP),U(IOUT),XO(IOUT),YO(IOUT),
     +          XS(IOUT),YS(IOUT),XD(IOUT),YD(IOUT),XDD(IOUT),YDD(IOUT)
      DIMENSION XP(IDIM),YP(IDIM),S(IDIM)
C
      DATA X/3.D0,4.D0,9.D0,16.D0,21.D0,27.D0,34.D0,36.D0,34.D0,26.D0,
     +     18.D0/
      DATA Y/2.4D0,9.6D0,14.4D0,12.0D0,9.6D0,8.4D0,13.2D0,21.6D0,30.0D0,
     +     37.2D0,38.4D0/
C
C  Do KURV1DP set up.
C
      SIGMA = 1.D0
      ISF = 3
      CALL KURV1DP(IDIM,X,Y,SLOP1,SLOP2,ISF,XP,YP,TEMP,S,SIGMA,IERR)
      IF (IERR.NE.0) THEN
          PRINT *,'Error return from KURV1DP =',IERR
          STOP
      END IF
C
C  Get interpolated points using KURV2DP.
C
      TINC = 1.0D0/ (IOUT-1)
      DO 10 I = 1,IOUT
          U(I) = (I-1)*TINC
          CALL KURV2DP(U(I),XO(I),YO(I),IDIM,X,Y,XP,YP,S,SIGMA)
   10 CONTINUE
C
C  Get the derivatives.
C
      DO 20 I = 1,IOUT
          CALL KURVDDP(U(I),XS(I),YS(I),XD(I),YD(I),XDD(I),YDD(I),IDIM,
     +                 X,Y,XP,YP,S,SIGMA)
   20 CONTINUE
C
C  Draw plot.
C
      CALL DRWFT4(IDIM,X,Y,IOUT,XO,YO,U,XD,YD)
C
      STOP
      END
      SUBROUTINE DRWFT4(II,X,Y,IOUT,XO,YO,U,XD,YD)
C
C  Define error file, Fortran unit number, and workstation type,
C  and workstation ID.
C
      PARAMETER (IERRF=6,LUNIT=2,IWTYPE=1,IWKID=1)
      DOUBLE PRECISION X
      DOUBLE PRECISION Y
      DOUBLE PRECISION XO
      DOUBLE PRECISION YO
      DOUBLE PRECISION U
      DOUBLE PRECISION XD
      DOUBLE PRECISION YD
C
C  Open GKS, open and activate a workstation.
C
      CALL GOPKS(IERRF,ISZDM)
      CALL GOPWK(IWKID,LUNIT,IWTYPE)
      CALL GACWK(IWKID)
C
C  Define a color table.
C
      CALL GSCR(IWKID,0,1.0D0,1.0D0,1.0D0)
      CALL GSCR(IWKID,1,0.0D0,0.0D0,0.0D0)
      CALL GSCR(IWKID,2,1.0D0,0.0D0,0.0D0)
      CALL GSCR(IWKID,3,0.0D0,1.0D0,0.0D0)
      CALL GSCR(IWKID,4,0.0D0,0.0D0,1.0D0)
C
C  Draw markers at original points.
C
      CALL BKGFT4(0.D0,40.D0,0.D0,40.D0,0.15D0,0.85D0,
     +            'Demo for KURV1DP/KURV2DP',0.035D0,0.5D0,0.93D0,0)
      CALL GRIDAL(4,5,4,5,1,1,10,0.D0,0.D0)
      CALL GSMKSC(2.D0)
      CALL GSPMCI(4)
      CALL GPM(II,X,Y)
C
C  Draw the interpolated curve
C
      CALL CURVE(XO,YO,IOUT)
      CALL FRAME
C
C  Plot the first derivatives of X and Y with respect to the parametric
C  variable U.
C
      CALL SET(0.D0,1.D0,0.D0,1.D0,0.D0,1.D0,0.D0,1.D0,1)
      CALL PCSETI('FN',21)
      CALL PLCHHQ(0.5D0,0.95D0,'Derivatives from KURVDDP',0.035D0,0.D0,
     +            0.D0)
      CALL BKGFT4(0.D0,1.D0,-80.D0,80.D0,0.55D0,0.87D0,'dx/du',0.030D0,
     +            0.65D0,0.82D0,1)
      CALL GRIDAL(5,5,4,5,1,1,10,0.D0,-80.D0)
      CALL CURVE(U,XD,IOUT)
      CALL BKGFT4(0.D0,1.D0,-40.D0,80.D0,0.10D0,0.42D0,'dy/du',0.030D0,
     +            0.39D0,0.37D0,1)
      CALL GRIDAL(5,5,3,5,1,1,10,0.D0,-40.D0)
      CALL CURVE(U,YD,IOUT)
      CALL FRAME
C
      CALL GDAWK(IWKID)
      CALL GCLWK(IWKID)
      CALL GCLKS
C
      RETURN
      END
      SUBROUTINE BKGFT4(XL,XR,YB,YT,YPB,YPT,LABEL,SIZL,POSXL,POSYL,IZL)
      DOUBLE PRECISION XL
      DOUBLE PRECISION XR
      DOUBLE PRECISION YB
      DOUBLE PRECISION YT
      DOUBLE PRECISION YPB
      DOUBLE PRECISION YPT
      DOUBLE PRECISION SIZL
      DOUBLE PRECISION POSXL
      DOUBLE PRECISION POSYL
      DOUBLE PRECISION XX
      DOUBLE PRECISION YY
      DIMENSION XX(2),YY(2)
      CHARACTER*(*) LABEL
C
      CALL SET(0.D0,1.D0,0.D0,1.D0,0.D0,1.D0,0.D0,1.D0,1)
      CALL PCSETI('FN',21)
      CALL PLCHHQ(POSXL,POSYL,LABEL,SIZL,0.D0,0.D0)
      CALL SET(0.17D0,0.87D0,YPB,YPT,XL,XR,YB,YT,1)
      IF (IZL.NE.0) THEN
          XX(1) = XL
          XX(2) = XR
          YY(1) = 0.D0
          YY(2) = 0.D0
          CALL GSPLCI(2)
          CALL GPL(2,XX,YY)
          CALL GSPLCI(1)
      END IF
C
      CALL GASETI('LTY',1)
      CALL PCSETI('FN',21)
      CALL GASETR('XLS',0.02D0)
      CALL GASETC('XLF','(F4.1)')
      CALL GASETR('YLS',0.02D0)
      CALL GASETC('YLF','(F5.1)')
      CALL GASETR('XMJ',0.02D0)
      CALL GASETR('YMJ',0.02D0)
C
      RETURN
      END

home | contents | defs | params | procedures | exmpls